【POJ】2888 Magic Bracelet

题意：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。

#include<cstdio>
#include<cstring>
#include<cmath>
#include<vector>
#define P 9973
#define EPS 1e-8
#define MAXN 12
#define MAXM 32000
using namespace std;
int n, m;
bool p[MAXM];
vector<int> prime;
vector<int> factor;
vector<int> primefactor;
struct Matrix {
    int mat[MAXN][MAXN];
    void Zero() {
        memset(mat, 0, sizeof(mat));
    }
    void Unit() {
        int i;
        for (i = 0; i < MAXN; i++)
            mat[i][i] = 1;
    }
} g;
Matrix operator *(Matrix &a, Matrix &b) {
    Matrix tmp;
    int i, j, k;
    tmp.Zero();
    for (k = 0; k < m; k++) {
        for (i = 0; i < m; i++) {
            if (a.mat[i][k] == 0)
                continue;
            for (j = 0; j < m; j++) {
                tmp.mat[i][j] += a.mat[i][k] * b.mat[k][j] % P;
                if (tmp.mat[i][j] >= P)
                    tmp.mat[i][j] -= P;
            }
        }
    }
    return tmp;
}
Matrix operator ^(Matrix a, int k) {
    Matrix tmp;
    tmp.Zero();
    tmp.Unit();
    for (; k; k >>= 1) {
        if (k & 1)
            tmp = tmp * a;
        a = a * a;
    }
    return tmp;
}
void Init() {
    int i, j;
    prime.clear();
    memset(p, true, sizeof(p));
    for (i = 2; i < 180; i++) {
        if (p[i]) {
            for (j = i * i; j < MAXM; j += i)
                p[j] = false;
        }
    }
    for (i = 2; i < MAXM; i++) {
        if (p[i])
            prime.push_back(i);
    }
}
int ExtGcd(int a, int b, int &x, int &y) {
    int t, d;
    if (b == 0) {
        x = 1;
        y = 0;
        return a;
    }
    d = ExtGcd(b, a % b, x, y);
    t = x;
    x = y;
    y = t - a / b * y;
    return d;
}
int Invmod(int a, int n) {
    int x, y;
    ExtGcd(a, n, x, y);
    return (x % n + n) % n;
}
void Factor(int x) {
    int i, tmp;
    factor.clear();
    tmp = (int) (sqrt((double) x) + EPS);
    for (i = 1; i <= tmp; i++) {
        if (x % i == 0) {
            factor.push_back(i);
            if (i == tmp && i * i == x)
                continue;
            factor.push_back(n / i);
        }
    }
}
int NoChange(int x) {
    int i, ans;
    Matrix tmp;
    tmp = g ^ x;
    for (i = ans = 0; i < m; i++) {
        ans += tmp.mat[i][i];
        if (ans >= P)
            ans -= P;
    }
    return ans;
}
void Prime(int x) {
    int i, tmp;
    primefactor.clear();
    tmp = (int) (sqrt((double) x) + EPS);
    for (i = 0; prime[i] <= tmp; i++) {
        if (x % prime[i] == 0) {
            primefactor.push_back(prime[i]);
            while (x % prime[i] == 0)
                x /= prime[i];
        }
    }
    if (x > 1)
        primefactor.push_back(x);
}
int Mul(int x, int &k) {
    int i, ans;
    ans = 1;
    for (i = k = 0; x; x >>= 1, i++) {
        if (x & 1) {
            k++;
            ans *= primefactor[i];
        }
    }
    return ans;
}
int Count(int x) {
    int i, j, t, ans, tmp;
    Prime(x);
    ans = 0;
    t = (int) primefactor.size();
    for (i = 1; i < (1 << t); i++) {
        tmp = Mul(i, j);
        if (j & 1)
            ans += x / tmp;
        else
            ans -= x / tmp;
    }
    return (x - ans) % P;
}
int Burnside() {
    int i, ans;
    Factor(n);
    for (i = ans = 0; i < (int) factor.size(); i++) {
        ans += Count(n / factor[i]) * NoChange(factor[i]) % P;
        if (ans >= P)
            ans -= P;
    }
    return ans * Invmod(n, P) % P;
}
int main() {
    int c;
    int i, j, k, x, y;
    Init();
    scanf("%d", &c);
    while (c--) {
        scanf("%d%d%d", &n, &m, &k);
        g.Zero();
        for (i = 0; i < m; i++) {
            for (j = 0; j < m; j++)
                g.mat[i][j] = 1;
        }
        while (k--) {
            scanf("%d%d", &x, &y);
            x--, y--;
            g.mat[x][y] = g.mat[y][x] = 0;
        }
        printf("%d\n", Burnside());
    }
    return 0;
}