http://www.cnblogs.com/wenruo/p/5304698.html
先看 Polya定理,Burnside引理 回忆一下基础知识。总结的很棒。
一个置换就是集合到自身的一个双射,置换群就是元素为置换的群。
再看 Polya入门 涨涨姿势。
Burnside定理,在每一种置换群也就是等价群中的数量和除以置换群的数量,即非等价的着色数等于在置换群中的置换作用下保持不变的着色平均数。
Polya定理:设G={π1,π2,π3........πn}是X={a1,a2,a3.......an}上一个置换群,用m中颜色对X中的元素进行涂色,那么不同的涂色方案数为:1/|G|*(mC(π1)+mC(π2)+mC(π3)+...+mC(πk)). 其中C(πk)为置换πk的循环节的个数。
代码:
ll gcd(ll a, ll b) { return b == 0 ? a : gcd(b, a % b); } ll pow(ll x, ll n) { ll res = 1; while (n) { if (n & 1) res *= x; x = x * x; n >>= 1; } return res; } ll polya(ll m, ll n) { ll ans = 0; for (ll i = 0; i < n; ++i) { ans += pow(m, gcd(n, i)); } if (n & 1) ans += n * pow(m, n / 2 + 1); else ans += (pow(m, n / 2) + pow(m, n / 2 + 1)) * n / 2; return ans / 2 / n; }
优化后
ll pow(ll x, ll n) { ll res = 1; while (n) { if (n & 1) res = res * x % p; x = x * x % p; n >>= 1; } return res; } ll eular(ll n) { ll res = 1; for (ll i = 2; i * i <= n; ++i) { if (n % i == 0) { n /= i; res = res * (i - 1); while (n % i == 0) { n /= i; res = res * i; } } } if (n > 1) res = res * (n - 1); return res % p; } ll polya(int m, int n) { ll sum = 0; ll i; for (i = 1; i * i < n; ++i) { if (n % i == 0) { sum += eular(i) * pow(m, n / i) % p; sum += eular(n / i) * pow(m, i) % p; } } if (i * i == n) sum += eular(i) * pow(m, i) % p; return sum / n; }
再优化一下欧拉函数
#define N 100000 int prime[N]; bool is_prime[N]; int sieve(int n) { int p = 0; for (int i = 0; i <= n; ++i) is_prime[i] = true; is_prime[0] = is_prime[1] = false; for (int i = 2; i <= n; ++i) { if (is_prime[i]) { prime[p++] = i; for (int j = 2 * i; j <= n; j += i) is_prime[j] = false; } } return p; } int phi(int n) { int rea = n; for(int i = 0; prime[i] * prime[i] <= n; i++) { if(n % prime[i] == 0) { rea = rea - rea / prime[i]; while (n % prime[i] == 0) n /= prime[i]; } } if(n > 1) rea = rea - rea / n; return rea; } ll polya(int m, int n) { ll sum = 0; ll i; for (i = 1; i * i < n; ++i) { if (n % i == 0) { sum += phi(i) * pow(m, n / i); sum += phi(n / i) * pow(m, i); } } if (i * i == n) sum += phi(i) * pow(m, i); if (n & 1) sum += n * pow(m, n / 2 + 1); else sum += (pow(m, n / 2) + pow(m, n / 2 + 1)) * n / 2; return sum / 2 / n; }