题意：m种珠子串成由n个珠子构成的环，并且有一些限制，比如第i种与第j种不能相邻，旋转后相同算是同一种方案，求方案数。(m<=10,n<=10^9)
如果没有限制条件，可以很容易用Burnside定理+容斥得到答案。
如果只考虑限制，观察发现，m很小，n很大。那么可以通过快速幂得到从第i种到第j种，共有k个珠子的方案数。
再回到Burnside定理，显然有n种置换，现在问题是如何求每种置换能保持不变的着色方案数：
在【POJ】2154 Color已经推出，第i种置换（即旋转360/n*i度）会有GCD(n,i)种不同的颜色，循环节长度为GCD(n,i)。
就可以用快速幂得到从第i种回到第i种，共有GCD(n,i)个珠子的方案数。带回Burnside定理+容斥，问题就解决了。
【POJ】2409 Let it Bead->【POJ】2154 Color->【POJ】2888 Magic Bracelet。
1 #include<cstdio> 2 #include<cstring> 3 #include<cmath> 4 #include<vector> 5 #define P 9973 6 #define EPS 1e-8 7 #define MAXN 12 8 #define MAXM 32000 9 using namespace std; 10 int n, m; 11 bool p[MAXM]; 12 vector<int> prime; 13 vector<int> factor; 14 vector<int> primefactor; 15 struct Matrix { 16 int mat[MAXN][MAXN]; 17 void Zero() { 18 memset(mat, 0, sizeof(mat)); 19 } 20 void Unit() { 21 int i; 22 for (i = 0; i < MAXN; i++) 23 mat[i][i] = 1; 24 } 25 } g; 26 Matrix operator *(Matrix &a, Matrix &b) { 27 Matrix tmp; 28 int i, j, k; 29 tmp.Zero(); 30 for (k = 0; k < m; k++) { 31 for (i = 0; i < m; i++) { 32 if (a.mat[i][k] == 0) 33 continue; 34 for (j = 0; j < m; j++) { 35 tmp.mat[i][j] += a.mat[i][k] * b.mat[k][j] % P; 36 if (tmp.mat[i][j] >= P) 37 tmp.mat[i][j] -= P; 38 } 39 } 40 } 41 return tmp; 42 } 43 Matrix operator ^(Matrix a, int k) { 44 Matrix tmp; 45 tmp.Zero(); 46 tmp.Unit(); 47 for (; k; k >>= 1) { 48 if (k & 1) 49 tmp = tmp * a; 50 a = a * a; 51 } 52 return tmp; 53 } 54 void Init() { 55 int i, j; 56 prime.clear(); 57 memset(p, true, sizeof(p)); 58 for (i = 2; i < 180; i++) { 59 if (p[i]) { 60 for (j = i * i; j < MAXM; j += i) 61 p[j] = false; 62 } 63 } 64 for (i = 2; i < MAXM; i++) { 65 if (p[i]) 66 prime.push_back(i); 67 } 68 } 69 int ExtGcd(int a, int b, int &x, int &y) { 70 int t, d; 71 if (b == 0) { 72 x = 1; 73 y = 0; 74 return a; 75 } 76 d = ExtGcd(b, a % b, x, y); 77 t = x; 78 x = y; 79 y = t - a / b * y; 80 return d; 81 } 82 int Invmod(int a, int n) { 83 int x, y; 84 ExtGcd(a, n, x, y); 85 return (x % n + n) % n; 86 } 87 void Factor(int x) { 88 int i, tmp; 89 factor.clear(); 90 tmp = (int) (sqrt((double) x) + EPS); 91 for (i = 1; i <= tmp; i++) { 92 if (x % i == 0) { 93 factor.push_back(i); 94 if (i == tmp && i * i == x) 95 continue; 96 factor.push_back(n / i); 97 } 98 } 99 } 100 int NoChange(int x) { 101 int i, ans; 102 Matrix tmp; 103 tmp = g ^ x; 104 for (i = ans = 0; i < m; i++) { 105 ans += tmp.mat[i][i]; 106 if (ans >= P) 107 ans -= P; 108 } 109 return ans; 110 } 111 void Prime(int x) { 112 int i, tmp; 113 primefactor.clear(); 114 tmp = (int) (sqrt((double) x) + EPS); 115 for (i = 0; prime[i] <= tmp; i++) { 116 if (x % prime[i] == 0) { 117 primefactor.push_back(prime[i]); 118 while (x % prime[i] == 0) 119 x /= prime[i]; 120 } 121 } 122 if (x > 1) 123 primefactor.push_back(x); 124 } 125 int Mul(int x, int &k) { 126 int i, ans; 127 ans = 1; 128 for (i = k = 0; x; x >>= 1, i++) { 129 if (x & 1) { 130 k++; 131 ans *= primefactor[i]; 132 } 133 } 134 return ans; 135 } 136 int Count(int x) { 137 int i, j, t, ans, tmp; 138 Prime(x); 139 ans = 0; 140 t = (int) primefactor.size(); 141 for (i = 1; i < (1 << t); i++) { 142 tmp = Mul(i, j); 143 if (j & 1) 144 ans += x / tmp; 145 else 146 ans -= x / tmp; 147 } 148 return (x - ans) % P; 149 } 150 int Burnside() { 151 int i, ans; 152 Factor(n); 153 for (i = ans = 0; i < (int) factor.size(); i++) { 154 ans += Count(n / factor[i]) * NoChange(factor[i]) % P; 155 if (ans >= P) 156 ans -= P; 157 } 158 return ans * Invmod(n, P) % P; 159 } 160 int main() { 161 int c; 162 int i, j, k, x, y; 163 Init(); 164 scanf("%d", &c); 165 while (c--) { 166 scanf("%d%d%d", &n, &m, &k); 167 g.Zero(); 168 for (i = 0; i < m; i++) { 169 for (j = 0; j < m; j++) 170 g.mat[i][j] = 1; 171 } 172 while (k--) { 173 scanf("%d%d", &x, &y); 174 x--, y--; 175 g.mat[x][y] = g.mat[y][x] = 0; 176 } 177 printf("%d\n", Burnside()); 178 } 179 return 0; 180 }
