题意简述
有一个n个点的有向图，每个点有两条出边，称为A边和B边。称一种构图是好的，当且仅当对于所有i，从第i个点出发，先连走x次A边，走1次B边，再连走y次A边，走1次B边，再走z次A边，最后回到了i。对合法的构图计数，模1e9+7。两个构图不同当且仅当某个点的A边或B边的终点不同。
\(1 \le n \le 1000\)
解法
很好的置换题。
把A边和B边分开，分成两个图。先特判x=y=z=0的情况，因为这时候A边是可以随便连的，而B边构成的图必须是一些环，且每个环的大小都是1或2。和其他的是不一样的。
考虑其余的情况。不妨令n个人分别在原图的n个点上，按照题目中说的走法同时走。发现不允许任意两个人在中途走到同一个点上，因为这样他们就无法分开，也无法到达各自的目的地了。所以如果令一个点A边指向的点为\(a_i\)，B边指向的点为\(b_i\)，则a和b都是1~n的一个排列。要知道排列(也就是置换)是可以像多项式一样相乘的，比如a和b相乘等于c，则\(c[i]=b[a[i]]\)。置换的相乘是满足结合律的，但是不满足交换律。 其实这道题要求的就是满足\(a^x \cdot b \cdot a^y \cdot b \cdot a^z=I\)的(a,b)的个数，其中\(I\)为排列\(\{1,2,\cdots n\}\)。
- 定义：\(a^{-1}\)是置换a的逆，满足\(a^{-1}[a[i]]=i。a \cdot a^{-1} = a^{-1} \cdot a =I\)
化简上面式子：
\[\begin{align} a^x \cdot b \cdot a^y \cdot b \cdot a^z&=I\\ a^{-x} \cdot a^x \cdot b \cdot a^y \cdot b \cdot a^z \cdot a^{-z} &= a^{-x}a^{-z}\\ b \cdot a^y \cdot b&=a^{-x-z}\\ \end{align} \]令\(c=ba^y\)，则\(c^2=a^{y-x-z}\)。可以把统计(a,b)的数量改成统计(a,c)的数量，这两个是等价的。如果\(y-x-z<0\)，就把式子右边变成\((a^{-1})^{y-x-z}\)，并统计(\(a^{-1},c\))的数量，和上面也是等价的。
令k=\(|y-x-z|\)。综合来讲，现在要干的事情就是对满足\(c^2=d^k\)的排列(c,d)计数。考虑确定了c，\(c^2\)会是什么样。对于一个排列(置换)，把\(i \to p_i\)连边，会得到一些环。对于c中的每个环，如果环长是奇数，则\(c^2\)中它还是一个环；否则会分裂成两个长度相等的环。再考虑确定了d，\(d^k\)会是什么样。发现对于d中一个长度为len的环，它会分裂成\(gcd(k,len)\)个长度为\(\frac{len}{gcd(k,len)}\)的环(手动画一下可知)。
终于要分析怎么计数了。令\(c^2=d^k=e\)。考虑最后e长什么样。对e的样子进行dp，\(dp_{i,j}\)表示现在要尝试往e中加入长度为i的环，当前已经用掉了j个节点的方案数。转移的时候直接枚举长度为i的环放k个。这样枚举的复杂度是调和级数，也就是\(O(n^2logn)\)，是没有问题的。假设现在已经枚举好了i,j,k，考虑怎么计算转移系数。先在剩下的\(n-j\)个点中选出\(i\cdot k\)个点，这是个组合数。然后把这些选出的点分成无顺序的k组，每组i个，然后再把每组排成一个有向环。这些都是简单的排列组合，不解释了。
然后重头戏是，这k个长度为i的环分别有多少种方案排成\(c\)和\(d\)。其实在dp中枚举每个i的时候，对这两种方案数分别做背包dp就可以了，这两个dp统共的复杂度是\(\sum_{i=1}^n(\frac ni)^2\)，反正是不大的。这里以对若干长度为i的环，有多少种方案排成\(d\)的dp举例。我们需要对每个长度len，预处理 会分裂成若干长度为len的小环的 所有可能的 环长值(之前说过d取了k次方是会分裂成小环的)。先枚举长度为i的环的个数k。然后设\(f_p\)表示已经使用了i个环中的p个的方案数。转移时先枚举 当前还没用的顺序最靠前的长度为i的环 所在的大环的长度(之前预处理过的)，然后系数用排列组合算一下即可。
总时间复杂度\(O(n^2logn+\sum_{i=1}^n(\frac ni)^2)\)。据说不用置换，直接在图上做也是可以的，但我不会。
#include <bits/stdc++.h>
#define rep(i,n) for(int i=0;i<n;++i)
#define repn(i,n) for(int i=1;i<=n;++i)
#define LL long long
#define pii pair <int,int>
#define fi first
#define se second
#define mpr make_pair
#define pb push_back
void fileio()
{
#ifdef LGS
freopen("in.txt","r",stdin);
freopen("out.txt","w",stdout);
#endif
}
void termin()
{
#ifdef LGS
std::cout<<"\n\nPROGRAM TERMINATED";
#endif
exit(0);
}
using namespace std;
const LL MOD=1e9+7;
LL qpow(LL x,LL a)
{
LL res=x,ret=1;
while(a>0)
{
if(a&1) (ret*=res)%=MOD;
a>>=1;
(res*=res)%=MOD;
}
return ret;
}
LL n,x,y,z,fac[2010],inv[2010],k,dp[2010][2010],combdp[2010][2010],combdp2[2010],f[2010],g[2010],pwii[2010];
vector <LL> cantot[2010];
LL C(LL nn,LL mm){return fac[nn]*inv[mm]%MOD*inv[nn-mm]%MOD;}
int main()
{
fileio();
freopen("game.in","r",stdin);
freopen("game.out","w",stdout);
fac[0]=1;repn(i,2005) fac[i]=fac[i-1]*i%MOD;
rep(i,2003) inv[i]=qpow(fac[i],MOD-2);
cin>>n>>x>>y>>z;
if(x==0&&y==0&&z==0)
{
LL ans=qpow(n,n),res=0;
rep(two,n/2+1)
{
LL val=C(n,two*2),full=two*2;
for(LL i=full-1;i>0;i-=2) (val*=i)%=MOD;
(res+=val)%=MOD;
}
(ans*=res)%=MOD;
cout<<ans<<endl;
termin();
}
k=llabs(x+z-y);
repn(len,n) cantot[len/__gcd(k,(LL)len)].pb(__gcd(k,(LL)len));
dp[0][0]=1;
rep(i,n)
{
LL ii=i+1,chomx=(n+ii-1)/ii;
rep(j,chomx+3) rep(k,chomx+3) combdp[j][k]=0;
combdp[0][0]=1;
rep(j,chomx) rep(k,j+1) if(combdp[j][k]>0)
{
(combdp[j+1][k+1]+=combdp[j][k])%=MOD;
if(k>0) (combdp[j+1][k-1]+=combdp[j][k]*k%MOD*ii)%=MOD;
}
rep(j,n+3) f[j]=0;
rep(j,chomx+1)
{
f[j]=0;
if(ii%2==0) f[j]=combdp[j][0];
else rep(k,j+1) (f[j]+=combdp[j][k])%=MOD;
}
pwii[0]=1;repn(j,n+3) pwii[j]=pwii[j-1]*ii%MOD;
rep(j,n+3) g[j]=0;
repn(jj,chomx)
{
rep(j,jj+2) combdp2[j]=0;
combdp2[0]=1;
rep(j,jj)
{
rep(k,cantot[ii].size()) if(j+cantot[ii][k]<=jj)
{
LL use=cantot[ii][k],mul=C(jj-j-1,use-1)*fac[use-1]%MOD*pwii[use-1]%MOD;
(combdp2[j+use]+=combdp2[j]*mul)%=MOD;
}
}
g[jj]=combdp2[jj];
}
rep(j,n+1) if(dp[i][j]>0)
{
(dp[i+1][j]+=dp[i][j])%=MOD;
LL pwc=1;
for(LL cho=1;cho*ii+j<=n;++cho)
{
(pwc*=inv[ii])%=MOD;
LL mul=C(n-j,cho*ii)*fac[cho*ii]%MOD*pwc%MOD*inv[cho]%MOD;
(mul*=qpow(fac[ii-1],cho))%=MOD;
(mul*=f[cho])%=MOD;(mul*=g[cho])%=MOD;
(dp[i+1][j+cho*ii]+=dp[i][j]*mul)%=MOD;
}
}
}
cout<<dp[n][n]<<endl;
termin();
}