poj2888
题意:给定N(10^9)种珠子及颜色总数M,并给出一些限制条件,如 颜色i与j不能相邻,求总共有多少种方案数,结果Mod 9973
思路:假设如果原题没有颜色限制的话,那就是一题polya的模板题,再加上phi(欧拉函数)的优化即可(雷同poj2154)
但是因为有了所谓的限制,那就不能直接套,我么该怎么办呢?
这我们很容易想到先前讲的那个以burnside的思想,polya的置换数的方法求,那怎么做?
我们把每个循环节数计算出来,那么在当前下,处于同一个环上的必须颜色相同,那么同一个环就可以看成一个点,然后对这些点进行染色,满足相邻的限制即可,这样的方案 数就是当前置换下的方案数。举个例子:
比如:n = 6, M = 3 , 限制为1与2不相邻
那么在旋转2个位置情况下,可得到置换:(1 3 5)(2 4 6)
那么就等价于用3种颜色给2个位置涂色,并满足颜色1与颜色2不相邻
这样的话就可以用矩阵来解决方案数的问题,把每种颜色间连一条边,如果颜色不能涂在相邻处便为0,否者为1,
设这个矩阵为A,那么A^m上的aij便表示i通过M个位置最后一个位置为J的方案数,然后再用快速幂优化即可。
PS: 这里还要涉及到形如 a/b mod c的问题
处理方法要用到费马小定理,因为 b与c互质,则 b ^(c -1)= 1 (mod c)
所以 a /b mod c= a /b *(b ^ (c - 1)) mod c = a* (b ^ (c - 2)) mod c
这样就可以化除为乘。。甚是神奇
code:
1 #include <iostream> 2 #include <cstring> 3 #include <cmath> 4 #include <algorithm> 5 #include <cstdio> 6 #include <cstdlib> 7 #define MXN 9973; 8 using namespace std; 9 10 int n , m, T, k; 11 int b[50][10][10], prime[20000], bo[50001], v[10][10], u[10][10], tot = 0, ans; 12 13 void find_prime(){ 14 for (int i = 2; i <= 50000; ++i) 15 if (!bo[i]){ 16 for (int j = i + i; j <= 50000; j += i) 17 bo[j] = 1; 18 prime[++tot] = i; 19 } 20 } 21 22 void multi_m(int a[][10], int b[][10] , int c[][10] ){ 23 for (int i = 0; i < m; ++i) 24 for (int j = 0; j < m; ++j){ 25 c[i][j] = 0; 26 for (int k = 0; k < m; ++k) 27 c[i][j] =(c[i][j] + a[i][k]* b[k][j]) % MXN; 28 } 29 30 /* for (int i = 0; i < m; ++i){ 31 for (int j = 0; j < m; ++j) 32 printf("%d ", c[i][j]); 33 puts(""); 34 } 35 puts("");*/ 36 } 37 38 void init(){ 39 int x, y; 40 scanf("%d%d%d",&n, &m, &k); 41 for (int i = 0; i < m; ++i) 42 for (int j = 0; j < m; ++j) 43 b[1][i][j] = 1; 44 for (int i = 1; i <= k; ++i){ 45 scanf("%d%d", &x, &y); 46 --x; --y; 47 b[1][x][y] = b[1][y][x] = 0; 48 } 49 for (int i = 2; i <= 30; ++i) 50 multi_m(b[i-1], b[i-1], b[i]); 51 52 /* for (int k = 1; k <= 30; ++k){ 53 for (int i = 0; i < m; ++i){ 54 for (int j = 0; j < m; ++j) 55 printf("%d ", b[k][i][j]); 56 puts(""); 57 } 58 puts("\n"); 59 }*/ 60 } 61 62 int phi(int now){ 63 int rea = now; 64 for (int i = 1; prime[i] * prime[i] <= now; ++i) 65 if (now % prime[i] == 0){ 66 rea -= rea / prime[i]; 67 while (now % prime[i] == 0) now /= prime[i]; 68 } 69 if (now > 1) rea -= rea / now; 70 return rea % MXN; 71 } 72 73 int cal(int now){ 74 for (int i = 0; i < m; ++i) 75 for (int j = 0; j < m; ++j) 76 if (i == j) v[i][j] = 1; 77 else v[i][j] = 0; 78 int i = 1, ret = 0; 79 80 while (now){ 81 if (now & 1) { 82 83 memcpy(u, v, sizeof(v)); 84 //memset(v, 0, sizeof(v)); 85 multi_m(u, b[i], v); 86 } 87 now >>= 1; 88 ++i; 89 } 90 for (int i = 0; i < m; ++i) 91 ret = (ret + v[i][i]) % MXN; 92 return ret; 93 } 94 95 void work(int now){ 96 int l = phi(now); 97 ans = (ans + l * cal(n / now)) % MXN; 98 } 99 100 void solve(){ 101 ans = 0; 102 for (int i = 1; i * i <= n; ++i) 103 if (n % i == 0){ 104 work(i); 105 if (i*i != n) work(n / i); 106 } 107 n %= MXN; 108 int num = MXN; 109 for (int i = 1; i <= num - 2; ++i) ans = (ans * n) % MXN; 110 printf("%d\n", ans); 111 } 112 113 int main(){ 114 freopen("poj2888.in", "r", stdin); 115 freopen("poj2888.out", "w", stdout); 116 scanf("%d", &T); 117 find_prime(); 118 while (T--){ 119 init(); 120 solve(); 121 } 122 fclose(stdin); fclose(stdout); 123 }