POJ 2888-Magic Bracelet(polya計數定理)

NO IMAGE

Magic Bracelet
Time Limit: 2000MS Memory Limit: 131072K
Total Submissions: 6024 Accepted: 1922
Description

Ginny’s birthday is coming soon. Harry Potter is preparing a birthday present for his new girlfriend. The present is a magic bracelet which consists of n magic beads. The are m kinds of different magic beads. Each kind of beads has its unique characteristic. Stringing many beads together a beautiful circular magic bracelet will be made. As Harry Potter’s friend Hermione has pointed out, beads of certain pairs of kinds will interact with each other and explode, Harry Potter must be very careful to make sure that beads of these pairs are not stringed next to each other.

There infinite beads of each kind. How many different bracelets can Harry make if repetitions produced by rotation around the center of the bracelet are neglected? Find the answer taken modulo 9973.

Input

The first line of the input contains the number of test cases.

Each test cases starts with a line containing three integers n (1 ≤ n ≤ 109, gcd(n, 9973) = 1), m (1 ≤ m ≤ 10), k (1 ≤ k ≤ m(m − 1) ⁄ 2). The next k lines each contain two integers a and b (1 ≤ a, b ≤ m), indicating beads of kind a cannot be stringed to beads of kind b.

Output

Output the answer of each test case on a separate line.

Sample Input

4
3 2 0
3 2 1
1 2
3 2 2
1 1
1 2
3 2 3
1 1
1 2
2 2
Sample Output

4
2
1
0
題目:POJ2888
題意:把n中珠子串成一個長度為m的環,有k個限制,被限制的兩種珠子不能相鄰。
思路:裸的polya計數加尤拉函式優化這裡就不再贅述如不懂的可以百度或者
這裡主要是如何求限制條件下的不重複答案。
假設我們迴圈節的長度為m,那麼我們要算長度為m的合法方案數的數。
考慮m<=10,那麼我們用一個矩陣記錄任意兩種珠子之間是否限制。a[i][j]=1表示第i種珠子和第j種珠子之間沒有限制。如果將這個矩陣自乘一次,那麼a[i][i]則表示∑a[i][j](1<=j<=m),舉個例子,如果m==4.
那麼自承一次後,a[4][4]=a[4][1] a[4][2] a[4][3] a[4][4],即長度為2的合法的方案數。按照這個規律,我們只要把這個矩陣n次方,就可以得到長度為n的合法方案的個數。
所以這裡配合矩陣快速冪就可以很快的求出我們要的答案了。
Tips:由於矩陣乘法的時候要取多次模,就會TLE,所以取模哪裡要優化一下。
AC程式碼:

#include<stdio.h>
#include<string.h>
#include<algorithm>
#include<queue>
#include<math.h>
#define met(s,k) memset(s,k,sizeof s)
#define scan(a) scanf("%d",&a)
#define scanl(a) scanf("%lld",&a)
#define scann(a,b) scanf("%d%d",&a,&b)
#define scannl(a,b) scanf("%lld%lld",&a,&b)
#define scannn(a,b,c) scanf("%d%d%d",&a,&b,&c)
#define prin(a) printf("%d\n",a)
#define prinl(a) printf("%lld\n",a)
using namespace std;
typedef long long ll;
const int maxn=50000;
const ll mod=9973;
int isprime[maxn],prime[maxn],pcont,n,m,k;
struct Matrix//模仿大牛的炫酷的矩陣快速冪寫法
{
ll mat[11][11];
friend Matrix operator* (Matrix x,Matrix y)
{
Matrix c;
met(c.mat,0);
for(int i=0;i<m;i  )
{
for(int j=0;j<m;j  )
{
for(int k=0;k<m;k  )
{
c.mat[i][j]=(c.mat[i][j] x.mat[i][k]*y.mat[k][j]);
}
c.mat[i][j]%=mod;//因為矩陣的係數不會很大,所以在計算過程中不會溢位,這裡把%放在迴圈外面,才不會TLE
}
}
return c;
}
friend Matrix operator^ (Matrix x,int y)
{
Matrix c;
met(c.mat,0);
for(int i=0;i<m;i  )c.mat[i][i]=1;
while(y)
{
if(y&1)c=c*x;
x=x*x;
y>>=1;
}
return c;
}
}A;
void init()
{
for(int i=2; i<maxn; i  )
{
if(!isprime[i])
{
prime[pcont  ]=i;
for(int j=2*i; j<maxn; j =i)
isprime[j]=1;
}
}
}
int ksm(int x,int y)//最後對答案要取逆元,也要用到快速冪
{
int res=1;
x%=mod;
while(y)
{
if(y&1)res=(res*x)%mod;
x=(x*x)%mod;
y>>=1;
}
return res;
}
int euler(int x)
{
int res=x;
for(int i=0;i<pcont&&prime[i]*prime[i]<=x;i  )
{
if(x%prime[i]==0)
{
res=res/prime[i]*(prime[i]-1);
while(x%prime[i]==0)x/=prime[i];
}
}
if(x!=1)res=res/x*(x-1);
return res%mod;
}
ll cal(int x)
{
ll res=0;
Matrix c=A^x;
for(int i=0;i<m;i  )res=(res c.mat[i][i])%mod;
return res;
}
ll polya()//裸的polya帶尤拉函式優化模板
{
ll ans=0;
for(int i=1; i*i<=n; i  )
{
if(n%i==0)
{
ans=(ans euler(i)*cal(n/i))%mod;
if(n/i!=i)ans=(ans euler(n/i)*cal(i))%mod;
}
}
return ans;
}
int main()
{
init();
int t;
scan(t);
while(t--)
{
scannn(n,m,k);
for(int i=0;i<m;i  )
for(int j=0;j<m;j  )A.mat[i][j]=1;
for(int i=0;i<k;i  )
{
int x,y;
scann(x,y);
A.mat[x-1][y-1]=0;
A.mat[y-1][x-1]=0;
}
ll ans=polya();
ans=(ans*ksm(n,mod-2))%mod;//求逆元
prinl(ans);
}
return 0;
}