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 也是一道极角排序+扫描的题。