http://poj.org/problem?id=2888

題意：Harry Potter 要用m種顏色的珠子做一個長度為n的手鐲，手鐲首尾相接。其中某些顏色的珠子不兼容，不能放在一起。求Harry Potter能夠早多少種不同的手鐲（每種顏色的珠子都有無限多顆，旋轉之后能夠吻合的算同一種）。

解法：polya計數，sum = sigma ( Euler( n / i )*Gettr( i ) ) / n    { i | n }     主要是計算Gettr( i )的問題我們把m種顏色的關系畫成一個無向圖，而Gettr(i)就是長度為 i 的回路的個數把無向圖表示成鄰接矩陣 G[maxn][maxn]，Gettr(i)就是這個矩陣的 i 次冪的跡，也就是 Trace(G^i)，G^i 可以用快速冪來求，可以先把 G^1 G^2 G^4 ... G^(2^31) 先預處理出來，加速快速冪的過程。

  1#include <cstdio>
  2#include <complex>
  3#include <cstdlib>
  4#include <iostream>
  5#include <cstring>
  6#include <string>
  7#include <algorithm>
  8using namespace std;
  9
 10const int maxn = 11;
 11const int maxm = 32000;
 12const int P = 9973;
 13int vis[ maxm ], p[ maxm ];
 14int plen, flen;
 15int aa[ 32 ], bb[ 32 ];
 16int n, m, ans;
 17int M[ 32 ][ maxn ][ maxn ];
 18void prime( )
 19{
 20    int i, j, k;
 21    plen = 0;
 22    memset( vis, 0, sizeof( vis ) );
 23    for( i = 2, k = 4; i < maxm; ++i, k += i + i - 1 )
 24    {
 25        if( !vis[ i ] )
 26        {
 27            p[ plen++ ] = i;
 28            if( k < maxm ) for( j = k; j < maxm; j += i ) vis[ j ] = 1;
 29        }
 30    }
 31}
 32
 33int pow_mod( int a, int b, int c )
 34{
 35    int ans = 1, d = a % c;
 36    while( b )
 37    {
 38        if( b & 1 ) ans = ans * d % c;
 39        d = d * d % c;
 40        b >>= 1;
 41    }
 42    return ans;
 43}
 44
 45void matrix_mul( int a[ ][ maxn ], int b[ ][ maxn ], int p )
 46{
 47    int c[ maxn ][ maxn ];
 48    for( int i = 0; i < m; i++ )
 49    {
 50        for( int j = 0; j < m; j++ )
 51        {
 52            c[ i ][ j ] = 0;
 53            for( int k = 0; k < m; k++ )
 54                c[ i ][ j ] += a[ i ][ k ] * b[ k ][ j ];
 55        }
 56    }
 57    for( int i = 0; i < m; i++ )
 58        for( int j = 0; j < m; j++ )
 59            a[ i ][ j ] = c[ i ][ j ] % p;
 60}
 61
 62int matrix_mod( int a[ ][ maxn ], int b, int c )
 63{
 64    int ans[ maxn ][ maxn ];
 65    memset( ans, 0, sizeof( ans ) );
 66    for( int i = 0; i < m; i++ ) ans[ i ][ i ] = 1;
 67    int j = 0;
 68    while( b )
 69    {
 70        if( b & 1 ) matrix_mul( ans, M[ j ], c );
 71        b >>= 1;
 72        j++;
 73    }
 74    int sum = 0;
 75    for( int i = 0; i < m; i++ )
 76        sum = ( sum + ans[ i ][ i ] ) % c;
 77    return sum;
 78}
 79
 80void factor( int n )
 81{
 82    flen = 0;
 83    for( int i = 0; p[ i ] * p[ i ] <= n; i++ )
 84    {
 85        if( n % p[ i ] == 0 )
 86        {
 87            for( bb[ flen ] = 0; n % p[ i ] == 0; ++bb[ flen ], n /= p[ i ] );
 88            aa[ flen++ ] = p[ i ];
 89        }
 90    }
 91    if( n > 1 ) bb[ flen ] = 1, aa[ flen++ ] = n;
 92} 
 93
 94int elur( int n )
 95{
 96    int phi = n;
 97    for( int i = 0; p[ i ] * p[ i ] <= n; i++ )
 98    {
 99        if( n % p[ i ] == 0 )
100        {
101            while( n % p[ i ] == 0 ) n /= p[ i ];
102            phi -= phi / p[ i ];
103        }
104    }
105    if( n > 1 ) phi -= phi / n;
106    return phi;
107}
108
109int inv( int n, int p )
110{
111    return pow_mod( n, elur( p ) - 1, p );
112}
113
114int e[ maxn ][ maxn ], g[ maxn ][ maxn ];
115
116void init( int c )
117{
118    for( int i = 0; i < m; i++ )
119        for( int j = 0; j < m; j++ )
120            M[ 0 ][ i ][ j ] = g[ i ][ j ];
121    for( int i = 1; i < 32; i++ )
122    {
123        for( int j = 0; j < m; j++ )
124        {
125            for( int k = 0; k < m; k++ )
126            {
127                M[ i ][ j ][ k ] = 0;
128                for( int l = 0; l < m; l++ )
129                    M[ i ][ j ][ k ] += M[ i - 1 ][ j ][ l ] * M[ i - 1 ][ l ][ k ];
130                M[ i ][ j ][ k ] %= c;
131            }
132        }
133    }
134}
135
136void DFS( int dep, int sum )
137{
138    if( dep == flen )
139    {
140        for( int i = 0; i < m; i++ )
141            for( int j = 0; j < m; j++ )
142                e[ i ][ j ] = g[ i ][ j ];
143        ans = ( ans + elur( n / sum ) % P * matrix_mod( e, sum, P ) % P ) % P;
144        return;
145    }
146    for( int tmp = 1, i = 0; i <= bb[ dep ]; tmp *= aa[ dep ], i++ )
147        DFS( dep + 1, sum * tmp );
148}
149
150int main(int argc, char *argv[])
151{
152    int k, x, y, cas;
153    prime( );
154    scanf( "%d", &cas );
155    while( cas-- )
156    {
157        scanf( "%d %d %d", &n, &m, &k );
158        for( int i = 0; i < m; i++ )
159            for( int j = 0; j < m; j++ )
160                g[ i ][ j ] = 1;
161        for( int i = 0; i < k; i++ )
162        {
163            scanf( "%d %d", &x, &y );
164            g[ x - 1 ][ y - 1 ] = g[ y - 1 ][ x - 1 ] = 0;
165        }
166        init( P );
167        factor( n );
168        ans = 0;
169        DFS( 0, 1 );
170        ans = ans * inv( n, P ) % P;
171        printf( "%d\n", ans ); 
172    }
173    return 0;
174}
175