[Burnside] POJ 2888

题意

在有珠子不能相邻的情况下,有几种方案

思路

orz见别人的题解&代码注释

代码

// https://blog.csdn.net/maxwei_wzj/article/details/73024349
#include <algorithm>
#include <bits/stdc++.h>
#include <cstdio>
#include <cstring>
#include <iostream>

using namespace std;
typedef long long ll;

const int INF = 0x7f7f7f7f;
const int maxn = 55555 + 10;
const int MOD = 9973;

int n, m, k; //串n个珠子,有m个颜色,每个颜色的数量不限
struct matrix {
    int a[ 11 ][ 11 ];

    friend matrix operator* ( matrix x, matrix y ) {
        matrix ret;
        memset ( &ret, 0, sizeof ( ret ) ); // memset
        for ( int i = 1; i <= m; ++i )
            for ( int j = 1; j <= m; ++j ) {
                for ( int k = 1; k <= m; ++k )
                    ret.a[ i ][ j ] += x.a[ i ][ k ] * y.a[ k ][ j ];
                ret.a[ i ][ j ] %= MOD; //矩阵乘法取模
            }
        return ret;
    }
} M;

int Pow ( int a, int b ) {
    int ans = 1;
    a %= MOD;

    while ( b ) {
        if ( b & 1 ) {
            ans = ( ans * a ) % MOD;
            --b;
        }
        b >>= 1;
        a = a * a % MOD;
    }
    return ans;
}

matrix Pow ( matrix a, int b ) {
    matrix ans;
    memset ( &ans, 0, sizeof ( ans ) );
    for ( int i = 1; i <= m; ++i )
        ans.a[ i ][ i ] = 1;
    while ( b ) {
        if ( b & 1 ) {
            ans = ans * a;
            --b;
        }
        b >>= 1;
        a = a * a;
    }
    return ans;
}

int cnt = 0;
int pri[ maxn ];
bool ispri[ maxn ];
void isPrime () {
    memset ( ispri, true, sizeof ( ispri ) );

    for ( int i = 2; i < maxn; ++i ) {
        if ( ispri[ i ] == true )
            pri[ cnt++ ] = i;

        for ( int j = i + i; j < maxn; j += i )
            ispri[ j ] = false;
    }
}

int phi ( int n ) {
    int ans = n;
    //枚举素数
    for ( int i = 0; pri[ i ] * pri[ i ] <= n; ++i ) {
        // n * ( 1 - 1/p1 ) ( 1- 1/p2 )...
        if ( n % pri[ i ] == 0 ) {
            ans = ans - ans / pri[ i ];

            while ( n % pri[ i ] == 0 )
                n /= pri[ i ];
        }
    }

    if ( n > 1 )
        ans = ans - ans / n;
    return ans % MOD;
}
/*
 * 现在需要求长度为i的子链有多少中可行的情况,且要保证头尾相连可行
 *f(i,k)代表长度为i的链子,最后的颜色是k
 *递推关系 f(i) = M^(i-1)*f(1)
 *现在需要求长度为d+1的链子,保证开头和结尾的颜色都是第k个
 *f(d+1) = M^d * f(1)
 *f(1)= [ m * 1 ],有m种情况,分别是第k位为1,其余位为0,代表第一个珠子取第k个颜色有1种
 *对每个f(1) = [ k位 = 1 ],对应每个f(d+1)取第k位的值,带入公式得到f(d+1,k) = M^d.[ k ][ k ]
 */
int getval ( matrix a ) {
    int ans = 0;
    for ( int i = 1; i <= m; ++i )
        ans += a.a[ i ][ i ];
    ans = ans % MOD;
    return ans;
}

void printmatrx ( matrix p ) {

    for ( int i = 1; i <= m; ++i ) {
        for ( int j = 1; j <= m; ++j )
            printf ( "%d ", p.a[ i ][ j ] );
        printf ( "\n" );
    }
}

int main () {
#ifdef LOCAL
    freopen ( "in", "r", stdin );
    // freopen("out","w",stdout);
#endif
    isPrime ();
    int t;
    scanf ( "%d", &t );
    while ( t-- ) {
        scanf ( "%d%d%d", &n, &m, &k );
        for ( int i = 1; i <= m; ++i )
            for ( int j = 1; j <= m; ++j )
                M.a[ i ][ j ] = 1;
        for ( int i = 0, u, v; i < k; ++i ) { //输入u,v
            scanf ( "%d%d", &u, &v );
            M.a[ u ][ v ] = M.a[ v ][ u ] = 0;
        }

        /*
         * burnside: 1/|G| * \Sigma ( |C(f)| ),|C(f)|代表对一个f,满足f * c = c的c的个数
         * polya:特殊条件下,n个珠子m个颜色,|C(f)| = m^(#k),#k指循环节的个数
         * 现在要求在有颜色限制条件下的C(f)
         * 旋转下保持相同,旋转i下的循环节个数 = gcd ( n,i ) = 需要排列颜色的长度
         * n = 6, i = 2, [ 1 3 5 ][ 2 4 6 ], gcd( n,i ) = 2,
         * 只需要排列1和2位,3 4,5 6和1 2保持一致,就能保证旋转后不变
         * 所以只需要计算长度为gcd( n,i )的子链有多少种颜色组合
         * 而gcd( n,i )有重复居多,欧拉函数优化
         */
        int ans = 0;
        for ( int i = 1; i * i <= n; ++i ) {
            if ( n % i == 0 ) {
                matrix tmp = Pow ( M, n / i );
                int val = getval ( tmp );
                ans = ( ans + phi ( i ) * val ) % MOD;

                if ( i * i != n ) {
                    tmp = Pow ( M, i );
                    val = getval ( tmp );
                    ans = ( ans + phi ( n / i ) * val ) % MOD;
                }
            }
        }
        int inv = Pow ( n % MOD, MOD - 2 ) % MOD;
        printf ( "%d\n", ans * inv % MOD );
    }

    return 0;
}