解題區►解題區 - 其他題目

HDU 5097 - Page Rank

Problem

配合上一篇巨量資料的學習，利用馬可夫過程找到 Page Rank。

Sample Input

1
2
3
4
5
4
0111
0011
0001
0100

Sample Output

1
0.15 1.49 0.83 1.53

Solution

然而 Page Rank 計算主要分成兩種，一種是指派所有的 Page Rank 總合為 S，如果加起來不足 S，則把不足的部分平攤給所有節點，這麼一來計算誤差就可以被碾平，但這樣會會因為太多網頁，而導致單一 Rank 過小。

另外一種，則單純倚靠隨機的全局擴散，雖然會有計算誤差，持續做下去一樣會收斂。

為了加快速度，將公式的 $1 - \beta$ 提出，沒有必要將一個稀疏矩陣變成稠密矩陣，如果每一個矩陣元素都加上 $\frac{1 - \beta}{N}$ 會導致整個矩陣不存在非 0 的元素，而事實上 0 的元素是相當多的。這麼一來每次迭代的效率為 $O(E)$，而非 $O(N^{2})$。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
#include <stdio.h>
#include <math.h>
#include <vector>
#include <string.h>
using namespace std;
// Markov process, math
const int MAXN = 4096;
const int MAXIT = 9999;
class PageRank {
public:
    const double eps = 1e-10;
    vector<int> g[MAXN], invg[MAXN];
    vector<double> r;
    int N;
    double beta, S;
    void init(int n, double beta) {
        this->N = n;
        this->beta = beta;
        for (int i = 0; i < n; i++)
            g[i].clear(), invg[i].clear();
    }
    void addEdge(int u, int v) {
        g[u].push_back(v);
        invg[v].push_back(u);
    }
    bool isComplete(vector<double> &a, vector<double> &b) {
        double e = 0;
        for (int i = (int) a.size() - 1; i >= 0; i--) {
            e += (a[i] - b[i]) * (a[i] - b[i]);
        }
        return e < eps;
    }
    void compute() {
        vector<double> next_r(N, 0);
        r.resize(N, 1.0);
        
        for (int it = 0; it < MAXIT; it++) {
            for (int i = 0; i < N; i++) {
                double tmp = 0;
                for (int j = 0; j < invg[i].size(); j++) {
                    int x = invg[i][j];
                    tmp += r[x] / g[x].size();
                }
                next_r[i] = tmp * beta + (1.0 - beta);
            }
            if (isComplete(r, next_r)) {
                r = next_r;
                return ;
            }
            r = next_r;
        }
    }
} g;
char s[MAXN][MAXN];
int main() {
    const double beta = 0.85f;
    int N;
    while (scanf("%d", &N) == 1) {
        for (int i = 0; i < N; i++)
            scanf("%s", s[i]);
        
        g.init(N, beta);
        for (int i = 0; i < N; i++) {
            for (int j = 0; j < N; j++) {
                if (s[i][j] == '1') {
                    g.addEdge(i, j);
                }
            }
        }
        
        g.compute();
        
        for (int i = 0; i < N; i++) {
            printf("%.2lf%c", g.r[i], i == N - 1 ? '\n' : ' ');
        }
    }
    return 0;
}
Read More +
學校課程►巨量資料

Page Rank 迭代計算 C++

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

1
2
3
4
5
6
7
8
9
10
11
0.7 3
3
2 1 2
1 2
1 2
   
0.85 3
3
2 1 2
1 2
1 0

Sample Output

1
2
3
4
5
6
7
A 0.300
B 0.405
C 2.295
   
A 1.163
B 0.644
C 1.192

Solution

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

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
#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
*/
Read More +
解題區►解題區 - UVa

UVa 11214 - Guarding the Chessboard

Problem

用最少的皇后，覆蓋盤面中所有的 ‘X’。皇后之間的相互攻擊可以忽略不理。

Sample Input

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
8 8
XXXXXXXX
XXXXXXXX
XXXXXXXX
XXXXXXXX
XXXXXXXX
XXXXXXXX
XXXXXXXX
XXXXXXXX
8 8
X.......
.X......
..X.....
...X....
....X...
.....X..
......X.
.......X
0

Sample Output

1
2
Case 1: 5
Case 2: 1

Solution

直接套 DLX 最小重複覆蓋的模板，但是單純使用會狂 TLE。額外加上初始貪心，找到搜索的第一把剪枝條件，這一個貪心使用每次找盤面上能覆蓋最多 ‘X’ 的位置，進行放置皇后的工作。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <vector>
#include <algorithm>
#include <assert.h>
using namespace std;
#define MAXV 0x3f3f3f3f
#define MAXE 1048576
#define MAXC 1048576
#define MAXR 65536
class DLX {
public:
    struct DacingLinks {
        int left, right;
        int up, down;
        int ch, rh;
        int data; // extra info
    } DL[MAXE];
    int s[MAXC], o[MAXR], head, size, Ans, findflag;
    void Remove(int c) {
        static int i;
        for(i = DL[c].down; i != c; i = DL[i].down) {
            DL[DL[i].right].left = DL[i].left;
            DL[DL[i].left].right = DL[i].right;
            s[DL[i].ch]--;
        }
    }
    void Resume(int c) {
        static int i;
        for(i = DL[c].down; i != c; i = DL[i].down) {
            DL[DL[i].right].left = i;
            DL[DL[i].left].right = i;
            s[DL[i].ch]++;
        }
    }
    int used[MAXC] = {};
    int H() {
        static int c, ret, i, j, time = 0;
        for(c = DL[head].right, ++time, ret = 0; c != head; c = DL[c].right) {
            if(used[c] != time) {
                ret ++, used[c] = time;
                for(i = DL[c].down; i != c; i = DL[i].down)
                    for(j = DL[i].right; j != i; j = DL[j].right)
                        used[DL[j].ch] = time;
            }
        }
        return ret;
    }
    void DFS(int k) {
        if(k + H() >= Ans) return;
        if(DL[head].right == head) {
            Ans = min(Ans, k);
            return;
        }
        int t = MAXV, c = 0, i, j;
        for(i = DL[head].right; i != head; i = DL[i].right) {
            if(s[i] < t) {
                t = s[i], c = i;
            }
        }
        for(i = DL[c].down; i != c; i = DL[i].down) {
            o[k] = i;
            Remove(i);
            for(j = DL[i].right; j != i; j = DL[j].right) Remove(j);
            DFS(k+1);
            for(j = DL[i].left; j != i; j = DL[j].left) Resume(j);
            Resume(i);
            if (findflag) break;
        }
    }
    int new_node(int up, int down, int left, int right) {
        assert(size < MAXE);
        DL[size].up = up, DL[size].down = down;
        DL[size].left = left, DL[size].right = right;
        DL[up].down = DL[down].up = DL[left].right = DL[right].left = size;
        return size++;
    }
    void addrow(int n, int Row[], int data) {
        int a, r, row = -1, k;
        for(a = 0; a < n; a++) {
            r = Row[a];
            DL[size].ch = r, s[r]++;
            DL[size].data = data;
            if(row == -1) {
                row = new_node(DL[DL[r].ch].up, DL[r].ch, size, size);
                DL[row].rh = a;
            }else {
                k = new_node(DL[DL[r].ch].up, DL[r].ch, DL[row].left, row);
                DL[k].rh = a;
            }
        }
    }
    void init(int m) {
        size = 0;
        head = new_node(0, 0, 0, 0);
        int i;
        for(i = 1; i <= m; i++) {
            new_node(i, i, DL[head].left, head);
            DL[i].ch = i, s[i] = 0;
        }
    }
} dlx;
const int dx[8] = {0, 0, 1, 1, 1, -1, -1, -1};
const int dy[8] = {1, -1, 0, 1, -1, 0, 1, -1};
int greedy(int n, int m, char g[][16]) {
    int ret = 0;
    while (true) {
        int has = 0;
        for (int i = 0; i < n; i++) {
            for (int j = 0; j < m; j++) {
                if (g[i][j] == 'X') {
                    has = 1;
                }
            }
        }
        if (!has)   break;
        int mx = 0, bx = 0, by = 0;
        for (int i = 0; i < n; i++) {
            for (int j = 0; j < m; j++) {
                int cnt = 0;
                for (int k = 0; k < 8; k++) {
                    int x = i + dx[k], y = j + dy[k];
                    while (x >= 0 && x < n && y >= 0 && y < m) {
                        if (g[x][y] == 'X')
                            cnt++;
                        x += dx[k], y += dy[k];
                    }
                }
                if (g[i][j] == 'X')
                    cnt++;
                if (cnt > mx)
                    mx = cnt, bx = i, by = j;
            }
        }
        for (int k = 0; k < 8; k++) {
            int x = bx + dx[k], y = by + dy[k];
            while (x >= 0 && x < n && y >= 0 && y < m) {
                if (g[x][y] == 'X')
                    g[x][y] = '.';
                x += dx[k], y += dy[k];
            }
        }
        if (g[bx][by] == 'X')
            g[bx][by] = '.';
        ret++;
    }
    return ret;
}
int main() {
    int cases = 0;
    int n, m;
    char g[16][16];
    while (scanf("%d %d", &n, &m) == 2 && n) {
        int cover_col = 0, A[16][16] = {};
        for (int i = 0; i < n; i++)
            scanf("%s", g[i]);
        for (int i = 0; i < n; i++) {
            for (int j = 0; j < m; j++) {
                if (g[i][j] == 'X') {
                    A[i][j] = ++cover_col;
                }
            }
        }
        dlx.init(cover_col);
   
        int row[2028], rowSize = 0;
        for (int i = 0; i < n; i++) {
            for (int j = 0; j < m; j++) {
                rowSize = 0;
                for (int k = 0; k < 8; k++) {
                    int x = i + dx[k], y = j + dy[k];
                    while (x >= 0 && x < n && y >= 0 && y < m) {
                        if (g[x][y] == 'X')
                            row[rowSize++] = A[x][y];
                        x += dx[k], y += dy[k];
                    }
                }
                if (g[i][j] == 'X')
                    row[rowSize++] = A[i][j];
                sort(row, row + rowSize);
                if (rowSize) {
                    dlx.addrow(rowSize, row, i * m + j);
                }
            }
        }
        dlx.Ans = greedy(n, m, g);
        dlx.DFS(0);
        printf("Case %d: %d\n", ++cases, dlx.Ans);
    }
    return 0;
}
Read More +
解題區►解題區 - UVa

UVa 1683 - In case of failure

Problem

對平面上，每一個點找到最鄰近點的距離，輸出歐幾里得距離的平方即可。

Sample Input

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
2 
10 
17 41 
0 34 
24 19 
8 28 
14 12 
45 5 
27 31 
41 11 
42 45 
36 27 
15 
0 0 
1 2 
2 3 
3 2 
4 0 
8 4 
7 4 
6 3 
6 1 
8 0 
11 0 
12 2 
13 1 
14 2 
15 0

Sample Output

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
200 
100 
149 
100 
149 
52 
97 
52 
360 
97 
5 
2 
2 
2 
5 
1 
1 
2 
4 
5 
5 
2 
2 
2 
5

Solution

KD Tree

計算幾何中的資料結構 KD-tree，雖然不算完美，一開始先進行 random search，先找到基礎的剪枝條件，隨後掛上啟發式進行搜索。

KD-tree 的效能在 D 度空間，區間詢問 (回報區域中的所有點或求區域極值) 的效能約為 O(n^(1-1/d))。我不會估算找鄰近點對的效能。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <algorithm>
#include <vector>
#include <queue>
using namespace std;
#define MAXN 262144
#define MAXD 2
#define INF (1LL<<60)
class kdTree {
public:
    struct PointD {
        long long d[MAXD];
        int label;
    };
    struct Node {
        Node *lson, *rson;
        int label;
    };
    struct cmp {
        bool operator()(const PointD* x, const PointD* y) const {
            if (x->d[sortDidx] != y->d[sortDidx])
                return x->d[sortDidx] < y->d[sortDidx];
            return x->label < y->label;
        }
    };
    struct pQcmp {
        bool operator() (pair<long long, int> &a, pair<long long, int> &b) const {
            if (a.first != b.first)	return a.first < b.first;
            return a.second < b.second;
        }
    };
    static int sortDidx;
    priority_queue<pair<long long, int>, vector< pair<long long, int> >, pQcmp> pQ; // <dist, label>
    Node buf[MAXN], *root;
    PointD pt[MAXN], *A[MAXN];
    int bufsize, K, qM, n;
    long long Q[MAXD], max_dist;
    void prepare(int n) {
        this->n = n;
        bufsize = 0;
        for (int i = 0; i < n; i++)
            A[i] = &pt[i];
        root = build(0, 0, n - 1);
    }
    Node* build(int k, int l, int r) {
        if(l > r)	return NULL;
        int m = (l + r)/2;
        Node *ret = &buf[bufsize++];
        sortDidx = k;
        nth_element(A + l, A + m, A + r + 1, cmp());
        ret->label = A[m]->label, ret->lson = ret->rson = NULL;
        if(l != r) {
            ret->lson = build((k+1)%K, l, m-1);
            ret->rson = build((k+1)%K, m+1, r);
        }
        return ret;
    }
    long long dist(int idx) {
        long long ret = 0;
        for(int i = 0; i < K; i++)
            ret += (long long)(pt[idx].d[i] - Q[i]) * (pt[idx].d[i] - Q[i]);
        return ret;
    }
    long long h_func(long long h[]) {
        long long ret = 0;
        for(int i = 0; i < K; i++)	ret += h[i];
        return ret;
    }
    void findNearest(Node *u, int k, long long h[]) {
        if(u == NULL || h_func(h) >= max_dist)
            return;
        long long d = dist(u->label);
        if (d <= max_dist) {
            pQ.push(make_pair(d, u->label));
            while (pQ.size() >= qM + 1)
                max_dist = pQ.top().first, pQ.pop();
        }
        long long old_hk = h[k];
        if (Q[k] <= pt[u->label].d[k]) {
            if (u->lson != NULL)
                findNearest(u->lson, (k+1)%K, h);
            if (u->rson != NULL) {
                h[k] = (pt[u->label].d[k] - Q[k]) * (pt[u->label].d[k] - Q[k]);
                findNearest(u->rson, (k+1)%K, h);
                h[k] = old_hk;
            }
        } else {
            if(u->rson != NULL)
                findNearest(u->rson, (k+1)%K, h);
            if (u->lson != NULL) {
                h[k] = (pt[u->label].d[k] - Q[k]) * (pt[u->label].d[k] - Q[k]);
                findNearest(u->lson, (k+1)%K, h);
                h[k] = old_hk;
            }
        }
    }
    void randomSearch(int ban) {
        max_dist = INF;
        long long d;
        for (int i = 0; i < 100; i++) {
            int x = rand()%n;
            if (x == ban)	continue;
            d = dist(x);
            if (d <= max_dist) {
                pQ.push(make_pair(d, x));
                while (pQ.size() >= qM + 1)
                    max_dist = pQ.top().first, pQ.pop();
            }
        }
    }
    void query(long long p[], int M, int ban) {
        for (int i = 0; i < K; i++)
            Q[i] = p[i];
        max_dist = INF, qM = M;
        long long h[MAXD] = {};
        randomSearch(ban);
        findNearest(root, 0, h);
        vector<int> ret;
        while (!pQ.empty())
            ret.push_back(pQ.top().second), pQ.pop();
        for (int i = (int) ret.size() - 1; i >= 0; i--) {
            if (ret[i] == ban)	continue;
            printf("%lld\n", dist(ret[i]));
            break;
        }
    }
};
int kdTree::sortDidx;
kdTree tree;
int main() {
    int n;
    long long p[MAXD];
    int testcase;
    scanf("%d", &testcase);
    while (testcase--) {
        scanf("%d", &n);
        tree.K = 2;
        for (int i = 0; i < n; i++) {
            for (int j = 0; j < 2; j++)
                scanf("%lld", &tree.pt[i].d[j]);
            tree.pt[i].label = i;
        }
        tree.prepare(n);
        for (int i = 0; i < n; i++) {
            for (int j = 0; j < 2; j++)
                p[j] = tree.pt[i].d[j];
            tree.query(p, 2, i);
        }
    }
    return 0;
}
/*
 2
 10
 17 41
 0 34
 24 19
 8 28
 14 12
 45 5
 27 31
 41 11
 42 45
 36 27
 15
 0 0
 1 2
 2 3
 3 2
 4 0
 8 4
 7 4
 6 3
 6 1
 8 0
 11 0
 12 2 
 13 1 
 14 2 
 15 0
 */

Delaunay

關於 Delaunay triangulation 在此就不說明，由於是 Voronoi Diagram，就可以找到所有相鄰的點。

效能是 O(n log n)，而平面上的邊為線性 O(n)，但這種做法必須離線。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
#include <stdio.h>
#include <math.h>
#include <algorithm>
#include <set>
#include <map>
#include <queue>
#include <vector>
#include <string>
#include <iostream>
#include <assert.h>
#include <string.h>
#include <list>
using namespace std;
#define eps 1e-8
#define MAXN (1048576)
#define MAXV 262144
struct Point {
    double x, y;
    int id;
    Point(double a = 0, double b = 0, int c = -1):
    x(a), y(b), id(c) {}
    Point operator-(const Point &a) const {
        return Point(x - a.x, y - a.y);
    }
    Point operator+(const Point &a) const {
        return Point(x + a.x, y + a.y);
    }
    Point operator*(const double a) const {
        return Point(x * a, y * a);
    }
    Point operator/(const double a) const {
        return Point(x / a, y / a);
    }
    bool operator<(const Point &a) const {
        if (fabs(x - a.x) > eps)	return x < a.x;
        if (fabs(y - a.y) > eps)	return y < a.y;
        return false;
    }
    bool operator==(const Point &a) const {
        return fabs(x - a.x) < eps && fabs(y - a.y) < eps;
    }
    bool operator!=(const Point &a) const {
        return !(fabs(x - a.x) < eps && fabs(y - a.y) < eps);
    }
    void read(int id = -1) {
        this->id = id;
    }
    double dist(Point b) {
        return hypot(x - b.x, y - b.y);
    }
    double dist2(Point b) {
        return (x - b.x) * (x - b.x) + (y - b.y) * (y - b.y);
    }
    void print() {
        printf("point (%lf, %lf)\n", x, y);
    }
};
struct Point3D {
    double x, y, z;
    Point3D(double a = 0, double b = 0, double c = 0):
    x(a), y(b), z(c) {}
    Point3D(Point p) {
        x = p.x, y = p.y, z = p.x * p.x + p.y * p.y;
    }
    Point3D operator-(const Point3D &a) const {
        return Point3D(x - a.x, y - a.y, z - a.z);
    }
    double dot(Point3D a) {
        return x * a.x + y * a.y + z * a.z;
    }
};
struct Edge {
    int id;
    list<Edge>::iterator twin;
    Edge(int id = 0) {
        this->id = id;
    }
};
int cmpZero(double v) {
    if (fabs(v) > eps)	return v > 0 ? 1 : -1;
    return 0;
}
double cross(Point o, Point a, Point b) {
    return (a.x-o.x)*(b.y-o.y) - (a.y-o.y)*(b.x-o.x);
}
Point3D cross(Point3D a, Point3D b) { // normal vector
    return Point3D(a.y * b.z - a.z * b.y
                   , -a.x * b.z + a.z * b.x
                   , a.x * b.y - a.y * b.x);
}
int inCircle(Point a, Point b, Point c, Point p) {
    if (cross(a, b, c) < 0)
        swap(b, c);
    Point3D a3(a), b3(b), c3(c), p3(p);
    //	printf("%lf %lf %lf\n", a3.x, a3.y, a3.z);
    //	printf("%lf %lf %lf\n", b3.x, b3.y, b3.z);
    //	printf("%lf %lf %lf\n", c3.x, c3.y, c3.z);
    //	printf("%lf %lf %lf\n", p3.x, p3.y, p3.z);
    b3 = b3 - a3, c3 = c3 - a3, p3 = p3 - a3;
    Point3D f = cross(b3, c3); // normal vector
    return cmpZero(p3.dot(f)); // check same direction, in: < 0, on: = 0, out: > 0
}
int intersection(Point a, Point b, Point c, Point d) { // seg(a, b) and seg(c, d)
    return cmpZero(cross(a, c, b)) * cmpZero(cross(a, b, d)) > 0
    && cmpZero(cross(c, a, d)) * cmpZero(cross(c, d, b)) > 0;
}
class Delaunay {
public:
    list<Edge> head[MAXV]; // graph
    Point p[MAXV];
    int n, rename[MAXV];
    void init(int n, Point p[]) {
        for (int i = 0; i < n; i++)
            head[i].clear();
        for (int i = 0; i < n; i++)
            this->p[i] = p[i];
        sort(this->p, this->p + n);
        for (int i = 0; i < n; i++)
            rename[p[i].id] = i;
        this->n = n;
        divide(0, n - 1);
    }
    void addEdge(int u, int v) {
        head[u].push_front(Edge(v));
        head[v].push_front(Edge(u));
        head[u].begin()->twin = head[v].begin();
        head[v].begin()->twin = head[u].begin();
    }
    void divide(int l, int r) {
        if (r - l <= 2) { // #point <= 3
            for (int i = l; i <= r; i++)
                for (int j = i+1; j <= r; j++)
                    addEdge(i, j);
            return;
        }
        int mid = (l + r) /2;
        divide(l, mid);
        divide(mid + 1, r);
        
        list<Edge>::iterator it;
        int nowl = l, nowr = r;
        
        //		printf("divide %d %d\n", l, r);
        for (int update = 1; update; ) { // find left and right convex, lower common tangent
            update = 0;
            Point ptL = p[nowl], ptR = p[nowr];
            for (it = head[nowl].begin(); it != head[nowl].end(); it++) {
                Point t = p[it->id];
                double v = cross(ptR, ptL, t);
                if (cmpZero(v) > 0 || (cmpZero(v) == 0 && ptR.dist2(t) < ptR.dist2(ptL))) {
                    nowl = it->id, update = 1;
                    break;
                }
            }
            if (update)	continue;
            for (it = head[nowr].begin(); it != head[nowr].end(); it++) {
                Point t = p[it->id];
                double v = cross(ptL, ptR, t);
                if (cmpZero(v) < 0 || (cmpZero(v) == 0 && ptL.dist2(t) < ptL.dist2(ptR))) {
                    nowr = it->id, update = 1;
                    break;
                }
            }
        }
        
        addEdge(nowl, nowr); // add tangent
        //		printf("add base %d %d\n", nowl, nowr);
        for (int update = 1; true;) {
            update = 0;
            Point ptL = p[nowl], ptR = p[nowr];
            int ch = -1, side = 0;
            for (it = head[nowl].begin(); it != head[nowl].end(); it++) {
                //				ptL.print(), ptR.print(), p[it->id].print();
                if (cmpZero(cross(ptL, ptR, p[it->id])) > 0
                    && (ch == -1 || inCircle(ptL, ptR, p[ch], p[it->id]) < 0))
                    ch = it->id, side = -1;
                //				printf("test L %d %d %d\n", nowl, it->id, inCircle(ptL, ptR, p[ch], p[it->id]));
            }
            for (it = head[nowr].begin(); it != head[nowr].end(); it++) {
                if (cmpZero(cross(ptR, p[it->id], ptL)) > 0
                    && (ch == -1 || inCircle(ptL, ptR, p[ch], p[it->id]) < 0))
                    ch = it->id, side = 1;
                //				printf("test R %d %d %d\n", nowr, it->id, inCircle(ptL, ptR, p[ch], p[it->id]));
            }
            if (ch == -1)	break; // upper common tangent
            //			printf("choose %d %d\n", ch, side);
            if (side == -1) {
                for (it = head[nowl].begin(); it != head[nowl].end(); ) {
                    if (intersection(ptL, p[it->id], ptR, p[ch])) {
                        head[it->id].erase(it->twin);
                        head[nowl].erase(it++);
                    } else
                        it++;
                }
                nowl = ch;
                addEdge(nowl, nowr);
            } else {
                for (it = head[nowr].begin(); it != head[nowr].end(); ) {
                    if (intersection(ptR, p[it->id], ptL, p[ch])) {
                        head[it->id].erase(it->twin);
                        head[nowr].erase(it++);
                    } else
                        it++;
                }
                nowr = ch;
                addEdge(nowl, nowr);
            }
        }
    }
    
    vector< pair<int, int> > getEdge() {
        vector< pair<int, int> > ret;
        list<Edge>::iterator it;
        for (int i = 0; i < n; i++) {
            for (it = head[i].begin(); it != head[i].end(); it++) {
                if (it->id < i)
                    continue;
                //				printf("DG %d %d\n", i, it->id);
                ret.push_back(make_pair(p[i].id, p[it->id].id));
            }
        }
        return ret;
    }
} tool;
#define INF (1LL<<60)
Point p[MAXV];
long long ret[MAXV];
int main() {
    int testcase, n;
    long long x, y;
    scanf("%d", &testcase);
    while (testcase--) {
        scanf("%d", &n);
        for (int i = 0; i < n; i++) {
            scanf("%lld %lld", &x, &y);
            p[i] = Point(x, y, i);
        }
        tool.init(n, p);
        vector< pair<int, int> > DG = tool.getEdge();
        for (int i = 0; i < n; i++) {
            ret[i] = INF;
        }
        for (int i = 0; i < DG.size(); i++) {
            x = DG[i].first, y = DG[i].second;
            long long v = (long long) (p[x].x - p[y].x) * (long long) (p[x].x - p[y].x) +
            (long long) (p[x].y - p[y].y) * (long long) (p[x].y - p[y].y);
            ret[x] = min(ret[x], v);
            ret[y] = min(ret[y], v);
        }
        for (int i = 0; i < n; i++) {
            printf("%lld\n", ret[i]);
        }
    }
    return 0;
}
/*
 2
 10
 17 41
 0 34
 24 19
 8 28
 14 12
 45 5
 27 31
 41 11
 42 45
 36 27
 15
 0 0
 1 2
 2 3
 3 2
 4 0
 8 4
 7 4
 6 3
 6 1
 8 0
 11 0
 12 2
 13 1
 14 2
 15 0
 */
Read More +
解題區►解題區 - UVa

UVa 1679 - Easy Geometry

Problem

在凸多邊形內部找一最大空白矩形 (平行兩軸)。多邊形頂點數最多 10 萬個。

Sample Input

1
2
3
4
5
6
7
8
9
10
11
4
5 1
2 4
3 7
7 3
5
1 1
1 4
4 7
7 4
7 1

Sample Output

1
2
3.5 2.5 5.5 4.5
1 1 7 4

Solution

三分內嵌三分再內嵌二分，幾何算法從函數觀點下手的有趣題目。

首先，觀察矩形在 x 軸的情況，觀察矩形左側 left.x，最大空白矩形面積呈現凸性，而在 left.x 下，right.x 分布也呈現凸性。

因此，三分左端點再三分右端點，最後二分查找兩條平行線交凸邊形的兩點，找到空白矩形的面積。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
#include <stdio.h> 
#include <math.h>
#include <algorithm>
#include <set>
#include <assert.h>
#include <vector>
using namespace std;
#define eps 1e-8
struct Pt {
    double x, y;
    Pt(double a = 0, double b = 0):
    	x(a), y(b) {}	
	Pt operator-(const Pt &a) const {
        return Pt(x - a.x, y - a.y);
    }
    Pt operator+(const Pt &a) const {
        return Pt(x + a.x, y + a.y);
    }
    Pt operator*(const double a) const {
        return Pt(x * a, y * a);
    }
    bool operator==(const Pt &a) const {
    	return fabs(x - a.x) < eps && fabs(y - a.y) < eps;
	}
    bool operator<(const Pt &a) const {
		if (fabs(x - a.x) > eps)
			return x < a.x;
		if (fabs(y - a.y) > eps)
			return y < a.y;
		return false;
	}
	double length() {
		return hypot(x, y);
	}
	void read() {
		scanf("%lf %lf", &x, &y);
	}
};
double dot(Pt a, Pt b) {
	return a.x * b.x + a.y * b.y;
}
double cross(Pt o, Pt a, Pt b) {
    return (a.x-o.x)*(b.y-o.y)-(a.y-o.y)*(b.x-o.x);
}
double cross2(Pt a, Pt b) {
    return a.x * b.y - a.y * b.x;
}
int between(Pt a, Pt b, Pt c) {
	return dot(c - a, b - a) >= -eps && dot(c - b, a - b) >= -eps;
}
int onSeg(Pt a, Pt b, Pt c) {
	return between(a, b, c) && fabs(cross(a, b, c)) < eps;
}
struct Seg {
	Pt s, e;
	double angle;
	int label;
	Seg(Pt a = Pt(), Pt b = Pt(), int l=0):s(a), e(b), label(l) {
//		angle = atan2(e.y - s.y, e.x - s.x);
	}
	bool operator<(const Seg &other) const {
		if (fabs(angle - other.angle) > eps)
			return angle > other.angle;
		if (cross(other.s, other.e, s) > -eps)
			return true;
		return false;
	}
	bool operator!=(const Seg &other) const {
		return !((s == other.s && e == other.e) || (e == other.s && s == other.e));
	}
};
Pt getIntersect(Seg a, Seg b) {
	Pt u = a.s - b.s;
    double t = cross2(b.e - b.s, u)/cross2(a.e - a.s, b.e - b.s);
    return a.s + (a.e - a.s) * t;
}
double getAngle(Pt va, Pt vb) { // segment, not vector
	return acos(dot(va, vb) / va.length() / vb.length());
}
Pt rotateRadian(Pt a, double radian) {
	double x, y;
	x = a.x * cos(radian) - a.y * sin(radian);
	y = a.x * sin(radian) + a.y * cos(radian);
	return Pt(x, y);
}
const double pi = acos(-1);
int cmpZero(double v) {
	if (fabs(v) > eps) return v > 0 ? 1 : -1;
	return 0;
}
#define MAXN 131072
#define INF 1e+30
Pt D[MAXN];
double ret;
Pt retLpt, retRpt;
double getArea(vector<Pt> &upper, vector<Pt> &lower, double lx, double rx) {
	int p, q;
	Pt a, b, c, d;
	p = (int) (lower_bound(upper.begin(), upper.end(), Pt(lx, INF)) - upper.begin());
	q = (int) (lower_bound(lower.begin(), lower.end(), Pt(lx, INF)) - lower.begin());
	assert(p > 0 && q > 0);
	a = getIntersect(Seg(Pt(lx, 1), Pt(lx, -1)), Seg(upper[p], upper[p-1]));
	b = getIntersect(Seg(Pt(lx, 1), Pt(lx, -1)), Seg(lower[q], lower[q-1]));
	p = (int) (lower_bound(upper.begin(), upper.end(), Pt(rx, INF)) - upper.begin());
	q = (int) (lower_bound(lower.begin(), lower.end(), Pt(rx, INF)) - lower.begin());
	assert(p > 0 && q > 0);
	c = getIntersect(Seg(Pt(rx, 1), Pt(rx, -1)), Seg(upper[p], upper[p-1]));
	d = getIntersect(Seg(Pt(rx, 1), Pt(rx, -1)), Seg(lower[q], lower[q-1]));
	double area = (rx - lx) * (min(a.y, c.y) - max(b.y, d.y));
	if (area > ret) {
		ret = area;
		retLpt = Pt(lx, max(b.y, d.y));
		retRpt = Pt(rx, min(a.y, c.y));
	}
	return area;
}
double search(vector<Pt> &upper, vector<Pt> &lower, double lx) {
	double l = lx, r = upper[upper.size()-1].x, mid, midmid, md, mmd;
	double ret = 0;
	for (int it = 0; it < 100; it++) {
		mid = (l + r)/2;
		midmid = (mid + r)/2;
		md = getArea(upper, lower, lx, mid);
		mmd = getArea(upper, lower, lx, midmid);
		ret = max(ret, md);
		if (md < mmd)
			l = mid;
		else
			r = midmid;
	}
	return ret;
}
double findMaxRectangle(vector<Pt> upper, vector<Pt> lower) {
	double mnx, mxx;
	mnx = upper[0].x, mxx = upper[upper.size()-1].x;
	double l = mnx, r = mxx, mid, midmid, md, mmd;
	double ret = 0;
	for (int it = 0; it < 100; it++) {
		mid = (l + r)/2;
		midmid = (mid + r)/2;
		md = search(upper, lower, mid);
		mmd = search(upper, lower, midmid);
		ret = max(ret, md);
		if (md < mmd)
			l = mid;
		else
			r = midmid;
	}
	return ret;
}
int main() {
	int N;
	while (scanf("%d", &N) == 1) {
		for (int i = 0; i < N; i++)
			D[i].read(); // clockwise order
		double mnx = INF, mxx = -INF;
		for (int i = 0; i < N; i++) {
			mnx = min(mnx, D[i].x);
			mxx = max(mxx, D[i].x);
		}
		vector<Pt> upper, lower;
		Pt Lpt = Pt(mnx, -INF), Rpt = Pt(mxx, INF);
		int lst, rst;
		for (int i = 0; i < N; i++) {
			if (cmpZero(D[i].x - mnx) == 0) {
				if (Lpt < D[i]) {
					Lpt = D[i], lst = i;
				}
			}
			if (cmpZero(D[i].x - mxx) == 0) {
				if (D[i] < Rpt) {
					Rpt = D[i], rst = i;
				}
			}
		}
		while (true) {
			upper.push_back(D[lst]);
			if (cmpZero(D[lst].x - mxx) == 0)
				break;
			lst = (lst + 1)% N;
		}
		while (true) {
			lower.push_back(D[rst]);
			if (cmpZero(D[rst].x - mnx) == 0)
				break;
			rst = (rst + 1)% N;
		}
		reverse(lower.begin(), lower.end());
		ret = 0;
		double area = findMaxRectangle(upper, lower);
//		printf("%lf\n", area);
		printf("%lf %lf %lf %lf\n", retLpt.x, retLpt.y, retRpt.x, retRpt.y);
	}
	return 0;
}
/*
4
5 1
2 4
3 7
7 3
5
1 1
1 4
4 7
7 4
7 1
*/
Read More +
解題區►解題區 - UVa

UVa 1676 - GRE Words Revenge

Problem

操作有二種，學習單詞以及給一篇文章，詢問學習單詞出現總次數。

輸入有經過加密，其加密的位移量為上一次詢問答案。

Sample Input

1
2
3
4
5
6
7
8
9
2
3
+01
+01
?01001
3
+01
?010
?011

Sample Output

1
2
3
4
5
Case #1:
2
Case #2:
1
0

Solution

如果是靜態，則可以想到 AC 自動機自然可以解決，但是 AC 自動機卻不支持動態建造。

為了滿足操作，可以利用快取記憶體的方式，進行分層，將小台 AC 自動機不斷重新建造，當小台 AC 自動機的節點總數大於某個值時，將其與大台 AC 自動機合併。合併操作必須是線性的，不要去重新插入每一個單詞。

特別小心有重複的單詞出現的小台或者是大台的操作，如果有重複必須忽略。雖然在理論上，分兩層操作的總複雜度高於倍增分層 (2, 4, 6, 8, …)，但是實作後，倍增分層一直 TLE。

然而這一題還有一種更困難的作法，詳見劉汝佳那本書，估計是一個動態後綴樹自動機 …

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
#include <iostream>
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <queue>
#include <map>
#define MAXCHAR 2
#define MAXS (5000005)
#define MAXNODE 526288
#pragma comment( linker, "/STACK:1024000000,1024000000")
using namespace std;
class ACmachine {
public:
    struct Node {
        Node *next[MAXCHAR], *fail;
        int cnt, val;
        void init() {
            fail = NULL;
            cnt = val = 0;
            memset(next, 0, sizeof(next));
        }
    } nodes[MAXNODE];
    Node *root;
    int size;
    Node* getNode() {
        Node *p = &nodes[size++];
        p->init();
        return p;
    }
    void init() {
        size = 0;
        root = getNode();
    }
    int toIndex(char c) {
        return c - '0';
    }
    void dfs(Node *p, Node *q) {
        for (int i = 0; i < MAXCHAR; i++) {
            if (q->next[i]) {
                Node *u = p->next[i], *v = q->next[i];
                if (u == NULL)
                    p->next[i] = getNode(), u = p->next[i];
                u->cnt |= v->cnt;
                dfs(u, v);
            }
        }
    }
    void merge(const ACmachine &other) {
        dfs(root, other.root);
    }
    void insert(const char str[]) {
        Node *p = root;
        for (int i = 0, idx; str[i]; i++) {
            idx = toIndex(str[i]);
            if (p->next[idx] == NULL)
                p->next[idx] = getNode();
            p = p->next[idx];
        }
        p->cnt = 1;
    }
    int find(const char str[]) {
        Node *p = root;
        for (int i = 0, idx; str[i]; i++) {
            idx = toIndex(str[i]);
            if (p->next[idx] == NULL)
                p->next[idx] = getNode();
            p = p->next[idx];
        }
        return p->cnt;
    }
    void build() { // AC automation
        queue<Node*> Q;
        Node *u, *p;
        Q.push(root), root->fail = NULL;
        while (!Q.empty()) {
            u = Q.front(), Q.pop();
            for (int i = 0; i < MAXCHAR; i++) {
                if (u->next[i] == NULL)
                    continue;
                Q.push(u->next[i]);
                p = u->fail;
                while (p != NULL && p->next[i] == NULL)
                    p = p->fail;
                if (p == NULL)
                    u->next[i]->fail = root;
                else
                    u->next[i]->fail = p->next[i];
                u->next[i]->val = u->next[i]->fail->val + u->next[i]->cnt;
            }
        }
    }
    int query(const char str[]) {
        Node *u = root, *p;
        int matched = 0;
        for (int i = 0, idx; str[i]; i++) {
            idx = toIndex(str[i]);
            while (u->next[idx] == NULL && u != root)
                u = u->fail;
            u = u->next[idx];
            u = (u == NULL) ? root : u;
            p = u;
            matched += p->val;
        }
        return matched;
    }
    void free() {
        return ;
        // owner memory pool version
        queue<Node*> Q;
        Q.push(root);
        Node *u;
        while (!Q.empty()) {
            u = Q.front(), Q.pop();
            for (int i = 0; i < MAXCHAR; i++) {
                if (u->next[i] != NULL) {
                    Q.push(u->next[i]);
                }
            }
            delete u;
        }
    }
};
char s[MAXS];
void decodeString(char s[], int L) {
    static char buf[MAXS];
    int n = (int) strlen(s);
    L %= n;
    for (int i = L, j = 0; i < n; i++, j++)
        buf[j] = s[i];
    for (int i = 0, j = n - L; i < L; i++, j++)
        buf[j] = s[i];
    for (int i = 0; i < n; i++)
        s[i] = buf[i];
}
ACmachine cache, disk;
int main() {
    int testcase, cases = 0;
    int N;
    scanf("%d", &testcase);
    while (testcase--) {
        scanf("%d", &N);
        printf("Case #%d:\n", ++cases);
        int L = 0; // encrypted key
        int cFlag = 0, dFlag = 0;
        cache.init(), disk.init();
        for (int i = 0; i < N; i++) {
            scanf("%s", s);
            decodeString(s+1, L);
            //            printf("%s\n", s);
#define THRESHOLD 4096
            if (s[0] == '+') {
                if (disk.find(s+1) || cache.find(s+1))
                    continue;
                cache.insert(s+1), cFlag = 0;
                if (cache.size > THRESHOLD) {
                    disk.merge(cache), dFlag = 0;
                    cache.free();
                    cache.init(), cFlag = 0;
                }
            } else {
                if (!cFlag) cache.build();
                if (!dFlag) disk.build();
                int ret = disk.query(s+1) + cache.query(s+1);
                printf("%d\n", ret);
                L = ret;
            }
        }
        disk.free(), cache.free();
    }
    return 0;
}
/*
 99999
 10
 +01
 +110
 ?010
 +110
 +00
 +0
 ?001001
 ?001001
 +110110
 ?1101001101
 
 6
 +01
 +110
 +110
 +00
 +0
 ?001001
 
 20
 +101001011
 ?110100
 +11010100
 ?0011001101
 +111011
 ?00010011
 +0111010110
 +0000101
 +0
 +11000
 ?1
 +1010101
 +0001
 +0110
 +0111101111
 ?1100
 +0111
 +1001
 ?0110111011
 ?1010010100
 
 10
 +00
 ?010110100
 +0100000100
 +111
 +000000
 ?0000110
 +110
 +00
 +0011
 ?101001
 
 5
 +0
 +1000100
 +01
 +0
 ?1110010011
 
 2
 3
 +01
 +01
 ?01001
 3
 +01
 ?010
 ?011
 */
Read More +
解題區►解題區 - UVa

UVa 1349 - Optimal Bus Route Design

Problem

每一條公車路線都是環狀，要求每一個站只在一條公車路線上，給定站與站之間的連通花費，求公車路線最少為何。

Sample Input

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
3 
2 2 3 1 0 
1 1 3 2 0 
1 3 2 7 0  
8 
2 3 3 1 0 
3 3 1 1 4 4 0
1 2 2 7 0 
5 4 6 7 0 
4 4 3 9 0 
7 4 8 5 0 
6 2 5 8 8 1 0
6 6 7 2 0
3 
2 1 0 
3 1 0 
2 1 0
0

Sample Output

1
2
3
7 
25 
N

Solution

看起來先拆點，要求每一個點變成出點跟入點，要求每一個出點都能配對一個入點，這麼一來能保證每一個點都能走出去，進而形成環狀。

那麼這就是一題完美帶權二分匹配判定。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
#include <stdio.h>
#include <string.h>
#include <queue>
#include <functional>
#include <deque>
#include <algorithm>
using namespace std;
#define MAXN 2048
#define MAXM 1048576
struct Node {
    int x, y, cap, cost;// x->y, v
    int next;
} edge[MAXM];
class MinCost {
public:
    int e, head[MAXN], dis[MAXN], pre[MAXN], record[MAXN], inq[MAXN];
    void Addedge(int x, int y, int cap, int cost) {
        edge[e].x = x, edge[e].y = y, edge[e].cap = cap, edge[e].cost = cost;
        edge[e].next = head[x], head[x] = e++;
        edge[e].x = y, edge[e].y = x, edge[e].cap = 0, edge[e].cost = -cost;
        edge[e].next = head[y], head[y] = e++;
    }
    pair<int, int> mincost(int s, int t) {
        int mncost = 0, flow, totflow = 0;
        int i, x, y;
        while(1) {
            memset(dis, 63, sizeof(dis));
            int oo = dis[0];
            dis[s] = 0;
            deque<int> Q;
            Q.push_front(s);
            while(!Q.empty()) {
                x = Q.front(), Q.pop_front();
                inq[x] = 0;
                for(i = head[x]; i != -1; i = edge[i].next) {
                    y = edge[i].y;
                    if(edge[i].cap > 0 && dis[y] > dis[x] + edge[i].cost) {
                        dis[y] = dis[x] + edge[i].cost;
                        pre[y] = x, record[y] = i;
                        if(inq[y] == 0) {
                            inq[y] = 1;
                            if(Q.size() && dis[Q.front()] > dis[y])
                                Q.push_front(y);
                            else
                                Q.push_back(y);
                        }
                    }
                }
            }
            if(dis[t] == oo)
                break;
            flow = oo;
            for(x = t; x != s; x = pre[x]) {
                int ri = record[x];
                flow = min(flow, edge[ri].cap);
            }
            for(x = t; x != s; x = pre[x]) {
                int ri = record[x];
                edge[ri].cap -= flow;
                edge[ri^1].cap += flow;
                edge[ri^1].cost = -edge[ri].cost;
            }
            totflow += flow;
            mncost += dis[t] * flow;
        }
        return make_pair(mncost, totflow);
    }
    void init(int n) {
        e = 0;
        for (int i = 0; i <= n; i++)
            head[i] = -1;
    }
} g;
int main() {
    int n, x, y, v;
    while (scanf("%d", &n) == 1 && n) {
        g.init(2 * n + 2);
        int source = 2 * n, sink = 2 * n + 1;
        for (int i = 0; i < n; i++) {
            x = i;
            while (scanf("%d", &y) == 1 && y) {
                scanf("%d", &v);
                y--;
                g.Addedge(2*x, 2*y+1, 1, v);
            }
            g.Addedge(source, 2*i, 1, 0);
            g.Addedge(2*i+1, sink, 1, 0);
        }
        pair<int, int> ret = g.mincost(source, sink);
        if (ret.second == n)
            printf("%d\n", ret.first);
        else
            puts("N");
    }
    return 0;
}
Read More +
解題區►解題區 - UVa

UVa 1343 - The Rotation Game

Problem

一款很神秘的遊戲，盤面上只會有 1, 2, 3 三種數字，可以藉由八種轉動方式，目標是將中間的八格變成都是相同數字。

求最少步數。

Sample Input

1
2
3
1 1 1 1 3 2 3 2 3 1 3 2 2 3 1 2 2 2 3 1 2 1 3 3
1 1 1 1 1 1 1 1 2 2 2 2 2 2 2 2 3 3 3 3 3 3 3 3
0

Sample Output

1
2
3
4
AC
2
DDHH
2

Solution

一開始設想，直接將三個分開討論，只考慮目標為 1 或 2 或 3，將狀態只變成 0/1，但分開計算的結果再取最小值可能會太慢，利用 BFS 很容易造成 TLE。

藉由 IDA* 的思路，啟發函數為中間位置還差多少個目標才能湊滿八個，因為每一次轉動，至多增加某一個數字在中間的個數。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
#include <stdio.h>
#include <map>
#include <set>
#include <queue>
#include <vector>
#include <string.h>
#include <sstream>
#include <iostream>
#include <algorithm>
using namespace std;
#define MAXSTATE 1048576
#define MAXN (20000000>>5)+1
#define GET(x) (mark[x>>5]>>(x&31)&1)
#define SET(x) (mark[x>>5] |= 1<<(x&31))
int mark[MAXN];
// 24! / (8! 16!) = 735471
int shift[8][8] = {
    {0, 2, 6, 11, 15, 20, 22},      // A
    {1, 3, 8, 12, 17, 21, 23},      // B
    {10, 9, 8, 7, 6, 5, 4},         // C
    {19, 18, 17, 16, 15, 14, 13},   // D
    {23, 21, 17, 12, 8, 3, 1},      // E
    {22, 20, 15, 11, 6, 2, 0},      // F
    {13, 14, 15, 16, 17, 18, 19},   // G
    {4, 5, 6, 7, 8, 9, 10}          // H
};
int center[8] = {6, 7, 8, 11, 12, 15, 16, 17};
int invShift[8] = {5, 4, 7, 6, 1, 0, 3, 2};
void shiftStrip(int A[], int dir) {
    int tmp = A[shift[dir][0]];
    for (int i = 0; i < 6; i++)
        A[shift[dir][i]] = A[shift[dir][i+1]];
    A[shift[dir][6]] = tmp;
}
int H(int A[]) {
    int K[4] = {};
    for (int i = 0; i < 8; i++) {
        K[A[center[i]]]++;
    }
    return 8 - max(K[1], max(K[2], K[3]));
}
int ida_depth, solved;
int path[128];
int IDA(int A[], int dep, int hv) {
    if (hv == 0) {
        solved = 1;
        if (dep == 0) {
            puts("No moves needed");
        } else {
            for (int i = 0; i < dep; i++)
                printf("%c", path[i] + 'A');
            puts("");
        }
        printf("%d\n", A[center[0]]);
        return dep;
    }
    if (dep + hv > ida_depth)
        return dep + hv;
    int back = 0x3f3f3f3f, shv, tmp;
    for (int i = 0; i < 8; i++) {
        shiftStrip(A, i);
        shv = H(A), path[dep] = i;
        tmp = IDA(A, dep + 1, shv);
        back = min(back, tmp);
        shiftStrip(A, invShift[i]);
        if (solved) return back;
    }
    return back;
}
int main() {
    while (true) {
        int A[32];
        for (int i = 0; i < 24; i++) {
            scanf("%d", &A[i]);
            if (A[i] == 0)  return 0;
        }
        solved = 0;
        for (ida_depth = 0; solved == 0;)
            ida_depth = IDA(A, 0, H(A));
    }
    return 0;
}
/*
 
*/
Read More +
解題區►解題區 - UVa

UVa 1336 - Fixing the Great Wall

Problem

有一台機器必須修復長城的所有位置，每一個位置的修復時間都是瞬間完成，但是移動必須耗時，而每一個位置的修復成本與修復時間呈現花費 c + d * t，求最少花費。

Sample Input

1
2
3
4
5
6
7
8
9
3 1 1000
1010 0 100
998 0 300
996 0 3
3 1 1000
1010 0 100
998 0 3
996 0 3
0 0 0

Sample Output

1
2
2084
1138

Solution

首先，可以藉由貪心得到，既然修復時間是瞬間，所經之處必然會修理完成。

因此狀態會呈現某一個區間已完成，並且停留在左端點或者是右端點，然而費用的計算則相當麻煩，因為計算某個區間而言，並不知道這區間花費的移動時間成本，只知道其最小花費。

反過來思考，利用潛熱的概念，預先支付往後所有未經過的位置所需要的花費，那麼就可以得到從左右端點彼此之間的移動時間 dt，將會增加未經過位置的所有 sumc，得到預支付花費 dt * sumc。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
#include <stdio.h> 
#include <math.h>
#include <string.h>
#include <algorithm>
using namespace std;
#define MAXN 1024
#define INF 1000000005
struct Pos {
	int x, c, delta; // cost = c + t * delta
	bool operator<(const Pos &a) const {
		return x < a.x;
	}
} D[MAXN];
int sumD[MAXN];
double dp[MAXN][MAXN][2]; // [interval_length][l][stop interval left/right]
int main() {
	int N, X;
	double V;
	while (scanf("%d %lf %d", &N, &V, &X) == 3 && N) {
		for (int i = 1; i <= N; i++)
			scanf("%d %d %d", &D[i].x, &D[i].c, &D[i].delta);
		sort(D + 1, D + 1 + N);
		
		sumD[0] = 0;
		for (int i = 1; i <= N; i++)
			sumD[i] = sumD[i-1] + D[i].delta;
		for (int i = 0; i <= N; i++) {
			for (int j = 0; j <= N + 1; j++) {
				dp[i][j][0] = dp[i][j][1] = INF;
			}
		}
		for (int i = 1; i <= N; i++) {
			double cost = sumD[N] * (fabs(X - D[i].x) / V) + D[i].c;
			dp[1][i][0] = dp[1][i][1] = cost;
		}
		
		for (int i = 2; i <= N; i++) {
			for (int j = 1; j + i - 1 <= N; j++) {
				int l = j, r = j + i - 1;
				double fromL, fromR;
				fromL = fromR = INF;
				fromL = dp[i-1][l+1][0] + (sumD[l] + sumD[N] - sumD[r]) * ((D[l+1].x - D[l].x) / V) + D[l].c;
				fromR = dp[i-1][l+1][1] + (sumD[l] + sumD[N] - sumD[r]) * ((D[r].x - D[l].x) / V) + D[l].c;
				dp[i][l][0] = min(fromL, fromR);
				
				fromL = dp[i-1][l][0] + (sumD[l-1] + sumD[N] - sumD[r-1]) * ((D[r].x - D[l].x) / V) + D[r].c;
				fromR = dp[i-1][l][1] + (sumD[l-1] + sumD[N] - sumD[r-1]) * ((D[r].x - D[r-1].x) / V) + D[r].c;
				dp[i][l][1] = min(fromL, fromR);
			}
		}
		printf("%d\n", (int) min(dp[N][1][0], dp[N][1][1]));
	}
	return 0;
}
Read More +
解題區►解題區 - UVa

UVa 1356 - Bridge

Problem

給吊橋總長、塔與塔之間的纜線限制、塔的高度，纜線在塔與塔之間呈現拋物線，使用最少的塔，請問纜線離地面多高。

Sample Input

1
2
3
4
5
2
20 101 400 4042
1 2 3 4

Sample Output

1
2
3
4
5
Case 1:
1.00
Case 2:
1.60

Solution

可以藉由貪心得到塔的數量，但是拋物線的方程式並不明確，因此可以考慮二分的高度，來計算纜線夠不夠長。

得到拋物線方程式後，要計算拋物線長度，使用數值積分的方式得到。為了降低計算誤差，使用自適應辛普森積分法，也就是當預估方法的變異小到一定程度時，停止分段求值。

integral parabola length

$$\int{\sqrt{dx^{2} + dy^{2}}} = \int{\sqrt{1 + dy^{2}/dx^{2}} dx} \\ = \int{\sqrt{1 + f'(x)^{2}} dx} $$ 
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
#include <stdio.h> 
#include <math.h>
using namespace std;
class SimpsonIntegral {
public:
	const double eps = 1e-6;
	double a; // function coefficient
	double f(double x) { // usually replace
		return sqrt(1 + 4 * a * a * x * x);
	}
	void setCoefficient(double a) {
		this->a = a;
	} 
	double simpson(double a, double b, double fa, double fb, double fab) {
		return (fa + 4 * fab + fb) * (b - a) / 6.0;
	}
	double integral(double l, double r) {
		return integral(l, r, eps);
	}
private:
	double integral(double a, double b, double c, double eps, double A, double fa, double fb, double fc) {
		double ab = (a + b)/2, bc = (b + c)/2;
		double fab = f(ab), fbc = f(bc);
		double L = simpson(a, b, fa, fb, fab), R = simpson(b, c, fb, fc, fbc);
		if (fabs(L + R - A) <= 15 * eps)
			return L + R + (L + R - A) / 15.0;
		return integral(a, ab, b, eps/2, L, fa, fab, fb) + integral(b, bc, c, eps/2, R, fb, fbc, fc);
	}
	double integral(double a, double b, double eps) {
		double ab = (a + b) /2;
		double fa = f(a), fb = f(b), fab = f(ab);
		return integral(a, ab, b, eps, simpson(a, b, fa, fb, fab), fa, fab, fb);
	}
} tool;
// integral parabola length
// \int \sqrt_{dx^{2} + dy^{2}} = \int \sqrt_{1 + dy^{2}/dx^{2}} dx
// = \int \sqrt_{1 + f'(x)^{2}} dx
int main() {
	int testcase, cases = 0;
	double D, H, B, L;
	scanf("%d", &testcase);
	while (testcase--) {
		scanf("%lf %lf %lf %lf", &D, &H, &B, &L);
		int n = ceil(B / D);
		double w = B / n, plen = L / n; // for a parabola
		double l = 0, r = H, len, mid;
		for (int it = 0; it < 100; it++) {
			mid = (l + r)/2; // parabola y = ax^2, pass (w/2, mid)
			tool.setCoefficient(4 * mid / w / w);
			len = 2 * tool.integral(0, w/2);
			if (len < plen)
				l = mid;
			else
				r = mid;
		}
		if (cases)	puts("");
		printf("Case %d:\n", ++cases);
		printf("%.2lf\n", H - l);
	}
	return 0;
}
Read More +