【JSWC2019】 小X的咒语

【JSWC2019】 小X的咒语

(\)

首先这道题有三个限制：

1. 每个点恰好两个出度和入度。
2. 没有自环。
3. 没有重边。

我们先定义几个变量：

(h_{i,j})：表示有(i)个出度入度为(1)的点和(j)个出度入度为(2)的点，可以有重边和自环的方案数。

(g_{i,j})：表示有(i+j)个出度入度为(2)的点，其中(i)个不能有自环，(j)个可以有自环的方案数。

(f_{i,j,k})：有(i)个点，这(i)个点之间只能连重边，不能连自环，连完边后还剩(j)个点没有连出边，其中的(k)个点也没有入边的方案数。

(ans_i)：至少有(i)个点连了重边，没有自环的方案数。

求出的(ans)过后就可以用容斥来算答案：

[Ans=sum_{i=0}^n(-1)^{i}ans_i ]

考虑求(ans_i)：

[ans_m=sum_{i=0}^{n-m}f_{n,n-m,i}*g_{i,n-m-i} ]

有(n-m)个点还没有连出边，假设这其中有(i)个点也没有入边，那么这(i)个点不能连出自环，另外的(n-m-i)的点不可能连出自环。

考虑求(f)：

用(DP)来求。考虑每次新加入一个点，它与之前的点的连边情况，具体方程就看代码吧。

考虑求(g)：

同样是用容斥。算至少有(i)个点连出自环的方案数。

[g_{i,j}=sum_{k=0}^i(-1)^{k}inom{i}{k}h_{k,i+j-k} ]

枚举(i)个点有自环，那么我们就先让这(i)个点都向自己连一条边，其他的乱连。

考虑求(h)：

假设有(n)个度数为(1)，(m)个度数为(2)的点。首先是((n+m*2)!)枚举每条边往哪里连，这样显然为算重。首先度数为(2)的点交换两条出边后仍是同种方案，所以先乘上(frac{1}{2^m})；交换两条入边也是相同的方案，所以在再乘上(frac{1}{2^m})就好了吗？不是的，因为如果某个点的两个入边是相同的点连过来的（也就是连出重边），那么就不能交换了。

设(G_i)表示至少有(i)个点连出重边，不考虑交换同一个点的两条入边的方案数，则：

[G_i=frac{1}{2^{m-i}}cdot inom{m}{i}cdot frac{m!}{(m-i)!}cdot (n+(m-i)*2)! ]

设(F_i)表示恰好(i)个点连出重边，不考虑交换同一个点的两条入边的方案数，则：

[F_i=sum_{j=i}(-1)^{j-i}inom{j}{i}G_j ]

然后：

[h=sum_{i=0}^mfrac{F_i}{2^{m-i}} ]

但是注意到，直接这么算是(O(N^4))的（枚举(n,m)，然后容斥算(F)又是(O(N^2))）。于是我们考虑将(h)写成：

[h=frac{sum_{i=0}^mcoef_icdot G_i}{2^m} ]

就是考虑每个(G)被计算的系数：

[coef_i=sum_{j=0}^i(-1)^{i-j}inom{i}{j}2^j ]

这个就是带入求(F)的容斥式子就可以推出来。

(\)

于是就做完了。

代码：

#include<bits/stdc++.h>
#define ll long long
#define N 505

using namespace std;
inline int Get() {int x=0,f=1;char ch=getchar();while(ch<'0'||ch>'9') {if(ch=='-') f=-1;ch=getchar();}while('0'<=ch&&ch<='9') {x=(x<<1)+(x<<3)+ch-'0';ch=getchar();}return x*f;}

int n,mod;
ll fac[N<<1],ifac[N<<1];
ll C[N][N],A[N][N];
ll pw[N<<1],pw2[N<<1];
ll inv2;
int cal(int i,int i1) {
	if(i1==0) return 1;
	else if(i1==1) return i;
	else return i*(i-1)/2;
}

ll sum[N];
ll coef[N];
void pre(int n) {
	fac[0]=1;
	for(int i=1;i<=2*n;i++) fac[i]=fac[i-1]*i%mod;
	for(int i=0;i<=n;i++)
		for(int j=0;j<=i;j++)
			C[i][j]=(!j||i==j)?1:(C[i-1][j-1]+C[i-1][j])%mod;
	for(int i=0;i<=n;i++) {
		A[i][0]=1;
		for(int j=1;j<=i;j++) A[i][j]=A[i][j-1]*(i-j+1)%mod;
	}
	inv2=mod+1>>1;
	pw[0]=pw2[0]=1;
	for(int i=1;i<=n;i++) {
		pw[i]=pw[i-1]*inv2%mod;
		pw2[i]=pw2[i-1]*2%mod;
	}
	for(int i=0;i<=n;i++) {
		for(int j=0;j<=i;j++) {
			if(i-j&1) {
				(coef[i]+=mod-pw2[j]*C[i][j]%mod+mod)%=mod;
			} else {
				(coef[i]+=pw2[j]*C[i][j]%mod)%=mod;
			}
		}
	}
}

//h: 乱连 
ll h[N][N];
//g:  无自环 
ll g[N][N];
//f:  
ll f[2][N][N];

ll cal_h(int n,int m) {
	static ll tem[N<<1],ans;
	for(int i=0;i<=m;i++) tem[i]=0;
	ans=0;
	for(int i=0;i<=m;i++) tem[i]=C[m][i]*A[m][i]%mod*fac[n+(m-i)*2]%mod*pw[m-i]%mod;

	for(int i=0;i<=m;i++) (ans+=tem[i]*coef[i])%=mod;
	return ans*pw[m]%mod;
}

ll Ans[N];
int main() {
	int type=Get();
	n=Get(),mod=Get();
	pre(n);
	if(n==1) {cout<<0;return 0;}
	for(int a=0;a<=n;a++)
		for(int b=0;b<=n-a;b++)
			h[a][b]=cal_h(a,b);
	for(int i=0;i<=n;i++) g[0][i]=h[0][i];
	for(int i=1;i<=n;i++) {
		for(int j=0;j<=n-i;j++) {
			g[i][j]=h[0][i+j];
			for(int k=1;k<=i;k++) {
				if(k&1) {
					g[i][j]=(g[i][j]-C[i][k]*h[k][i+j-k])%mod;
					if(g[i][j]<0) g[i][j]+=mod;
				} else {
					g[i][j]=(g[i][j]+C[i][k]*h[k][i+j-k])%mod;
				}
			}
			if(g[i][j]<0) g[i][j]+=mod;
			else if(g[i][j]>=mod) g[i][j]-=mod;
		}
	}
	f[0][0][0]=1;
	for(int i=0;i<n;i++) {
		int now=i&1;
		memset(f[now^1],0,sizeof(f[now^1]));
		for(int j=0;j<=i;j++) {
			for(int k=0;k<=j;k++) {
				if(!f[now][j][k]) continue ;
				(f[now^1][j+1][k+1]+=f[now][j][k]);
				int a=k,b=j-k,c=j-k;
				if(b) (f[now^1][j][k]+=f[now][j][k]*b);
				if(k) (f[now^1][j][k-1]+=f[now][j][k]*k);
				
				if(c) (f[now^1][j][k]+=f[now][j][k]*c);
				if(k) (f[now^1][j][k-1]+=f[now][j][k]*k);
				
				if(k) (f[now^1][j-1][k-1]+=f[now][j][k]*k);
				if(k>1) (f[now^1][j-1][k-2]+=f[now][j][k]*k*(k-1));
				if(k&&b) (f[now^1][j-1][k-1]+=f[now][j][k]*k*b);
				if(k&&c) (f[now^1][j-1][k-1]+=f[now][j][k]*k*c);
				if(b&&c) (f[now^1][j-1][k]+=f[now][j][k]*b*c);
			}
		}
		for(int j=0;j<=i+1;j++)
			for(int k=0;k<=j;k++)
				if(f[now^1][j][k]) f[now^1][j][k]%=mod;
	}
	int now=n&1;
	for(int i=0;i<=n;i++)
		for(int j=0;j<=i;j++) {
			(Ans[n-i]+=g[j][i-j]*f[now][i][j])%=mod;
		}
	for(int i=n;i>=0;i--)
		for(int j=i+1;j<=n;j++)
			Ans[i]=(Ans[i]-Ans[j]*C[j][i]%mod+mod)%mod;
	cout<<Ans[0];
	return 0;
}