Northwestern Europe 2005 Unequalled Consumption求最小的q, 使方程a1*x1+a2*x2+...+an*xn=q的整數解大于給定的數。
n<=5, ai <= 10,q <= 10^15
設該方程的解數是f(q)，則f(q)未必單調，但對確定的q0, f(q0+t*a1)必然關于t單調（因為t小的時候的所有的解都可以平行移動到大的t的解）
然后枚舉下q0（反正才10個），二分。
問題的關鍵就是求f(q)
利用母函數，f(q)就是母函數
P(x) = (1 + x^a1 + x^(2*a1) + x^(3*a1) +...)(1 + x^a2 + x^(2*a2) + ...)...(1 + x^an + x^(2*an) + ...)
的x^q項系數。
令LCM為a1, a2...an的最小公倍數。
則P(x) =
(1 + x^a1 + x^(2*a1) + ... + x^(LCM - a1)) * (1 + x^LCM + x^(2 * LCM) + ...)*
(1 + x^a2 + x^(2*a2) + ... + x^(LCM - a2)) * (1 + x^LCM + x^(2 * LCM) + ...)*
*...*
(1 + x^an + x^(2*an) + ... + x^(LCM - an)) * (1 + x^LCM + x^(2 * LCM) + ...)
=
(1 + x^a1 + x^(2*a1) + ... + x^(LCM - a1)) * (1 + x^a2 + x^(2*a2) + ... + x^(LCM - a2)) * ... *(1 + x^an + x^(2*an) + ... + x^(LCM - an)) *
(1 + x^LCM + x^(2 * LCM) + ...)^n
注意前一部分的x^r系數可以算出來（例如用DP）
后一項中x^(LCM * k)的系數是C(n+k-1, n-1) （推一下就知道了）
然后就可以得出結果了是吧～～～
p.s.
利用推出的結論可以知道對于給定的r, f(q + LCM * t)是一個關于t的n次函數。。。所以其實可以算出所有小于n * LCM的f(q)值然后使用Langrage插值公式，這個是標程的做法。（其實我想知道有沒有更簡單的方法證明這是一個關于t的n次函數～）
SPOJ 598, INCR求n的排列中，最長上升序列為b的排列個數。
n<=40, b <= 5
做法是dp, 狀態記錄n的排列中len = 2的上升序列最后元素的最早位置，len = 3的上升序列最后元素的最早位置...
每次轉移的時候枚舉n + 1的放置位置，并計算新的狀態。
由于有效狀態的稀疏性，可以用map / hash來優化。。。
CODE
1
// SPOJ598, INCR
2// Written By FreePeter
3
4#include <cstdio>
5#include <cstring>
6#include <iostream>
7#include <map>
8
9using namespace std;
10
11const int MaxN = 40 + 5, MOD = 1000000000;
12map<int, int> f[MaxN];
13
14void all_init();
15bool judge(int id, int b);
16int encode(int stat[4]);
17void decode(int id, int stat[4]);
18void make_newstat(int stat[4], int pos, int new_stat[4]);
19
20int main() 
{
21
all_init();
22
23 int t;
24 cin >> t;
25
for (; t > 0; --t) {
26 int ans = 0;
27 int n, b;
28 cin >> n >> b;
29 for (map<int, int>::iterator it = f[n].begin(); it != f[n].end(); ++it) {
30 if (judge(it->first, b)) {
31 ans += it->second;
32 if (ans >= MOD) ans -= MOD;
33
}
34 }
35
36 cout << ans << endl;
37 }
38
39 return 0;
40
}
41
42void all_init() {
43 int stat[4] = {50, 50, 50, 50};
44 f[1][encode(stat)] = 1;
45 for (int i = 1; i < 40; ++i) {
46 for (map<int, int>::iterator it = f[i].begin(); it != f[i].end(); ++it) {
47 decode(it->first, stat);
48
49 int new_stat[4];
50 for (int j = 0; j <= i; ++j) {
51 if (j > stat[3]) break; // It can make an permutation whose degree will exceed 5
52
53 make_newstat(stat, j, new_stat);
54 int &val = f[i + 1][encode(new_stat)];
55 val += it->second;
56 if (val >= MOD) val -= MOD;
57 }
58 }
59 }
60}
61
62bool judge(int id, int b) {
63 int stat[4];
64 decode(id, stat);
65 for (int i = 0; i < b - 1; ++i)
66 if (stat[i] == 50) return false;
67 for (int i = b - 1; i < 4; ++i)
68 if (stat[i] != 50) return false;
69
70 return true;
71}
72
73int encode(int stat[4]) {
74 int ans = 0;
75 for (int i = 0; i < 4; ++i)
76 ans = (ans << 6) | stat[i];
77 return ans;
78}
79
80void decode(int id, int stat[4]) {
81 for (int i = 3; i >= 0; --i) {
82 stat[i] = id & 63;
83 id >>= 6;
84 }
85}
86
87void make_newstat(int stat[4], int pos, int new_stat[4]) {
88 for (int i = 0; i < 4; ++i) {
89 if (stat[i] < pos) new_stat[i] = stat[i];
90 else {
91 if (i == 0) new_stat[i] = (0 < pos ? pos : min(stat[i] + 1, 50));
92 else new_stat[i] = (stat[i - 1] < pos ? pos : min(stat[i] + 1, 50));
93 }
94 }
95
96}
97
題意懶得說了～
首先要有1是吧。。。然后我們按照1的個數來分類，我們來計算恰有k個1的方案數。
我們在k個1的基礎上加入新的數，顯然第一個數只能是k+1
然后加入的數只能是k + 1 or 2 * (k + 1)
如法炮制。。。發現非1的數都具有(k + 1) * t的形式。。。設其依次為(k + 1) * ti
則{ti}這些數也滿足題目的性質。。。共有f((n - k) / (k + 1))種方案。
設f(n)是要求的函數，則f(n) = sigma(f((n - k) / (k + 1)), (k + 1) | (n + 1) , k>=1
f(0) = 1
這樣直接做會T...
我們令g(n) = f(n - 1)
則g(n) = f(n - 1) = sigma(f((n - k - 1) / (k + 1))) = sigma(f(n / (k + 1) - 1)) = sigma(g(n / (k + 1)), (k + 1) | n, k >= 1
設n = p1^a1 * p2^a2 * ... * pr*ar
令h(p1, p2,.., pr, a1, a2...ar) = g(n)
= h(p1,p2, ...pr, b1, b2, ...br),
0<=bi <= ai, bi不全=ai
注意對于一個確定的n，h()中的p1, p2...pr在計算過程中始終不變。。。所以。。。
計算結果與pi無關,只與ai有關這樣狀態數就大大減少了。。。直接因式分解后dp就行了。。。
CODE
1// Dhaka 2007 Problem C, The Dumb Grocer
2// Written By FreePeter
3
4#include <cstdio>
5#include <cstring>
6#include <iostream>
7#include <map>
8#include <algorithm>
9#include <limits>
10
11using namespace std;
12
13const int MaxN = 100000, MaxINT = numeric_limits<int>::max();
14map<long long, long long> f;
15long long fact[11];
16int prime[MaxN / 8], prime_cnt = 0;
17bool is_prime[MaxN + 10];
18
19long long encode(int stat[10]);
20void decode(long long id, int stat[10]);
21long long calc(long long id);
22long long generate_all(int deep, int stat[10], int cur_stat[10]);
23void make_prime();
24
25int main() {
26 make_prime();
27
28 fact[0] = 1;
29 for (int i = 1; i <= 10; ++i)
30 fact[i] = fact[i - 1] * i;
31
32 f[0] = 1;
33 int t;
34 cin >> t;
35 for (int testno = 0; testno < t; ++testno) {
36 int n;
37 cin >> n;
38 if (n == MaxINT) {
39 cout << "Case " << testno + 1 << ": " << n / 2 + 1 << endl;
40 continue;
41 }
42 ++n;
43
44 int e[10] = {0}, e_cnt = 0;
45 for (int i = 0; prime[i] * prime[i] <= n; ++i) {
46 if (n % prime[i] == 0) {
47 int exp = 1; n /= prime[i];
48 for (; n % prime[i] == 0; n /= prime[i]) ++exp;
49 e[e_cnt++] = exp;
50 }
51 }
52 if (n != 1) e[e_cnt++] = 1;
53
54 sort(e, e + 10);
55 cout << "Case " << testno + 1 << ": " << calc(encode(e)) << endl;
56 }
57
58 return 0;
59}
60
61long long encode(int stat[10]) {
62 long long ans = 0;
63 for (int i = 0; i < 10; ++i)
64 ans = (ans << 6) | stat[i];
65 return ans;
66}
67
68void decode(long long id, int stat[10]) {
69 for (int i = 9; i >= 0; --i) {
70 stat[i] = id & 63;
71 id >>= 6;
72 }
73}
74
75long long calc(long long id) {
76 map<long long, long long>::iterator it = f.find(id);
77 if (it != f.end()) return it->second;
78
79 int stat[10], cur_stat[10];
80 decode(id, stat);
81 long long ans = generate_all(0, stat, cur_stat);
82 return f[id] = ans;
83}
84
85long long generate_all(int deep, int stat[10], int cur_stat[10]) {
86 if (deep == 10) {
87 if (memcmp(cur_stat, stat, sizeof(int) * 10) == 0) return 0;
88
89 long long ans = 1;
90 int ptr = 9;
91 for (int i = 9; i >= 0; --i) {
92 for (; ptr >= 0 && stat[ptr] >= cur_stat[i]; ) --ptr;
93 int fac = 9 - ptr;
94 fac -= (9 - i);
95
96 ans *= fac;
97 }
98
99 ptr = 0;
100 for (; ptr < 10; ) {
101 int new_ptr = ptr;
102 for (; new_ptr < 10 && cur_stat[new_ptr] == cur_stat[ptr]; ) ++new_ptr;
103 ans /= fact[new_ptr - ptr];
104 ptr = new_ptr;
105 }
106
107 return ans *= calc(encode(cur_stat));
108 }
109
110 long long ans = 0;
111 for (int i = (deep == 0 ? 0 : cur_stat[deep - 1]); i <= stat[deep]; ++i) {
112 cur_stat[deep] = i;
113 ans += generate_all(deep + 1, stat, cur_stat);
114 }
115
116 return ans;
117}
118
119void make_prime() {
120 fill(is_prime, is_prime + MaxN + 10, true);
121 for (int i = 2; i <= MaxN; ++i)
122 if (is_prime[i]) {
123 prime[prime_cnt++] = i;
124 if (i >= 1000) continue;
125 for (int j = i * i; j <= MaxN; j += i)
126 is_prime[j] = false;
127 }
128}
129
130
Dhaka 2007 You are around me ...首先旋轉坐標，變成平行與xy軸的橢圓，然后坐標伸縮。。。變成圓。。。最近點對。。。貼模板。。。
ZJU2107 Quoit Design 一道測最近點對的題。
Dhaka 2007 Magnetic Train Tracks給定n個點，求可以構成多少個銳角三角形。
n <= 1200
話說求銳角三角形不太好算是吧。。。補集轉換，我們來求鈍角/直角三角形 <=> 求鈍角/直角個數。。。
后面的事情就簡單了，是對每個點，將其他點按照極角排序 + 掃描。
Dhaka 2005 Counting Triangles 也是一道補集轉換的題～（轉化成求三點共線的個數）
Shanghai 2004 Amphiphilic Carbon Molecules 也是一道極角排序+掃描的題。