Toy(生成函数+Burnside引理+矩阵快速幂)

On birthday, Anthony got a toy. It is constructed with N+1(N>=3) balls and 2*N sticks. All balls are in a same plane. One of them is special, while the other N balls are connected to it by N sticks with the same length. The angles between any two adjacent sticks are equal. And finally, any two adjacent balls(except the central one) are connected by a stick. 
  Here are two examples: 
 
Anthony wanted to remove N sticks, leaving all balls still connected. He wanted to know the number of all legal solutions. Your task is to solve this problem for him. 
  Notice that if a solution will be the same as another one by rotation, these two solutions should be consider as the same. 
The answer may be quite large. You just need to calculate the remainder of the answer when divided by M. 

Input

Input contains several test cases. 
For each test case, there is only one line containing two integers N and M(3<=N<=10^9, 2<=M<=10^9).
Input is terminated by EOF. 

Output

For each case, output one integer in one line, representing the remainder of the number of all solutions when divided by M. 
 

Sample Input

3 10000
4 10000
4 10

Sample Output

6
13
3

 

 

 

 

题解

难啊。。。

由于要考虑旋转等价,所以先上来一个Burnside引理

不动点的个数怎么算呢?

\sum_{i=1}^{n}F[n/gcd(n,i)]

先枚举一下置换,然后再算该置换下的不动点数

i表示当前置换把x置换到了x+i

那么置换的总个数就是gcd(n,i),每个置换的长度就是n/gcd(n,i)

设f[i]表示有i个点的轮状图的生成树方案数

那么当置换的总个数为gcd(n,i)时,我们可以把相邻的gcd(n,i)个点缩成一个点然后来求方案数

 

由于直接枚举的复杂度比较高,所以我们考虑枚举gcd

于是,不动点数就变成了:

\sum_{i|n}\varphi(n/i)*F[i](读者自行理解不难)

于是就可以在O(sqrt(n))的时间复杂度做了

然后让我们来考虑一下F[n]怎么计算

 

首先,我们假设没有环,就有一个比较明显的DP,记作f[n](我们把由外环边连通的点叫做一个块)

如图:(图中就有四个块,每个块都要派一个节点来与中心相连)

 

f[n]=1 (n=0、1)

f[n]=f[n-1](直接连向中心)+f[n-1](直接连向n-1点所在的块)+sumf(0~n-2)(把自己与中心相连,然后枚举前面连了多少个)

然后利用高考的混合型数列的化简技巧:

f[n]=f[n-1]+sum[n-1]

f[n-1]=f[n-2]+sum[n-2]

f[n]-f[n-1]=f[n-2]-f[n-2]+f[n-1](两式相减)

f[n]=3*f[n-1]-f[n-2]

 

然后我们再来考虑有环的情况

利用破环为链

我们可以枚举第一个点的所在块的大小,那么剩余的点就是链的情况了

设我们枚举1号点所在块的大小为i,则每个i都会贡献i个方案(因为一个大小为i的块指派节点的方案数有i种)

于是我们的最终答案F[n],则

F[n]=\sum_{i=1}^ni^2f[n-i]

长得很“卷积”

 

接下来就是鬼畜的生成函数推导时间

 

首先我们要弄明白f[n]的生成函数闭形式是什么

设f[n]=3*f[n-1]-f[n-2]的生成函数是f(x),那么f(x)就满足:

f(x)=3x*f(x)-x^2*f(x)

(1-3x+x^2)f(x)=0

!!!!怎么肥四?f(x) 的生成函数闭形式竟然是0

我们开始意识到问题不对劲。。。

原来f[n]=3*f[n-1]-f[n-2]递推式只有再n>=3的时候适用,f[0],f[1],f[2]并不满足

我们为了得到正确的式子生成函数不得不把首几项凑出来

f[0]=1,f[1]=1,f[2]=3,f[4]=8........(简写为1,1,3,8...)

f(x)   ->   1,1,3,8,21...

-3x*f(x)   ->   0,-3,-3,-9,-24...

+x^2*f(x)   ->   0,0,1,1,3...

f(x)-3x*f(x)+x^2*f(x)   ->   1,-2,1,0,0...

我们发现通过实际的f(x)来推出的等式右边并不等于0,但是第三项之后的系数都为0了

所以实际上f(x)的闭形式为:

f(x)=\frac{x^2-2x+1}{x^2-3x+1}

 

 

然后我们再来看g(x)=\sum_{i=0}^\infty i^2x^i的闭形式

我们已经知道了\sum_{i=0}^{\infty}x^i=\frac{1}{1-x},需要想办法凑出g(x)

原系数:1,1,1,1……

求导:1,2,3,4……

乘x:0,1,2,3……

求导:1,4,9,16……

再乘x:0,1,4,9,16……

成功了!!!( •̀ ω •́ )y

我们来看看它闭形式的变化:

原闭形式:\frac{1}{1-x}

求导:\frac{1}{(1-x)^2}

乘x:\frac{x}{(1-x)^2}

求导:\frac{1-x^2}{(1-x)^4}=\frac{1+x}{(1-x)^3}

再乘x:\frac{x+x^2}{(1-x)^3}

所以g(x)=\frac{x+x^2}{(1-x)^3}

 

 

激动人心的时刻到了:

F(x)=f(x)*g(x)=\frac{x^2-2x+1}{x^2-3x+1}*\frac{x+x^2}{(1-x)^3}

=\frac{x+x^2}{(1-x)*(x^2-3x+1)}

=\frac{x+x^2}{1-4x+4x^2-x^3}

移项:

F(x)-4x*F(x)+4x^2*F(x)-x^3F(x)=x+x^2

我们发现等式右边从x^3起都没有了系数

所以当n>=3时,F[n]=4*F[n-1]-4*F[n-2]+F[n-3]

然后把首项依次算出来:

F[0]=0,F[1]=1,F[2]=5

就可以愉快地矩阵快速幂啦!!!( •̀ ω •́ )y

 

 

肝了我三天才完全搞懂生成函数的推导。。。

代码:

#include<cstdio>
#include<cstring>
#include<algorithm>
using namespace std;
#define LL long long 
LL mo;
LL mul(LL x,LL y)
{
	LL t=0;x%=mo;y%=mo;
	while(y){
		if(y&1){t+=x;if(t>=mo)t-=mo;}
		y>>=1;x+=x;if(x>=mo)x-=mo;
	}
	return (t+mo)%mo;
}
struct mat{
	LL a[4][4];
	mat(){memset(a,0,sizeof(a));}
};
mat operator *(mat A,mat B)
{
	mat C;
	for(int i=1;i<=3;i++)
		for(int j=1;j<=3;j++)
			for(int k=1;k<=3;k++)
				C.a[i][j]=(C.a[i][j]+mul(A.a[i][k],B.a[k][j]))%mo;
	return C;
}
LL cal(int x)
{
	if(x==1) return 1;
	if(x==2) return 5;
	mat A,B;
	A.a[1][1]=5;A.a[2][1]=1;A.a[3][1]=0;
	B.a[1][1]=4;B.a[1][2]=mo-4;B.a[1][3]=1;
	B.a[2][1]=1;
	B.a[3][2]=1;
	x-=2;
	for(;x;x>>=1,B=B*B)
		if(x&1) A=B*A;
	return A.a[1][1];
}
int phi(int x)
{
	int ans=x;
	for(int i=2;i*i<=x;i++)if(x%i==0){
		ans=ans/i*(i-1);
		do{x/=i;}while(x%i==0);
	}
	if(x>1) ans=ans/x*(x-1);
	return ans%mo;
}
int main()
{
	LL ans;int n,m;
	while(~scanf("%d%d",&n,&m)){
		mo=1ll*n*m;ans=0;
		for(int i=1;i*i<=n;i++)if(n%i==0){
			ans=(ans+mul(phi(n/i),cal(i)))%mo;
			if(i*i!=n) ans=(ans+mul(phi(i),cal(n/i)))%mo;
		}
		printf("%d\n",(ans/n+m)%m);
	}
}

 

 

 

 

 

Burnside引理是最基本的,Polya定理只是Burnside引理在染色方案下计算的一个简便方法。

所以等价类计算的关键还是在于计算不动点

 

 

 

 

评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值