【gdgzezoi】Problem C: Shuffle and Swap

Description

你有两个长度都为N的01字符串A和B。A中的1的个数和B中的1的个数相等。

你想通过以下方式把A变成B：

1. 设a1,a2,…,ak是所有A中的1的下标。
2. 设b1,b2,…,bk是所有B中的1的下标。
3. 将a1,a2,…,ak随机排序，每种排列出现的概率都是1k!，这一步会产生k!种不同的结果。
4. 将b1,b2,…,bk随机排序，每种排列出现的概率都是1k!，这一步会产生k!种不同的结果。
5. 对于所有1≤i≤k，依次交换Aai和Abi。

记P为这样搞完之后A和B变得相等的概率。显然，P×(k!)2是个整数，你需要计算这个整数对998244353取模后的值。

Input
输入包含两行，每行一个01字符串，分别表示A和B。

Output
输出一个整数：P×(k!)2mod998244353。
Sample Input
Sample Input 1
1010
1100

Sample Input 2
01001
01001

Sample Input 3
101010
010101

Sample Input 4
1101011011110
0111101011101
Sample Output
Sample Output 1
3

Sample Output 2
4

Sample Output 3
36

Sample Output 4
932171449
HINT
1≤|A|=|B|≤104
A,B均由0、1组成，包含1的个数相等，且至少包含一个1。

存在40分的子任务，满足1≤|A|=|B|≤500

思路

我们可以把答案拆分成两步：

1.枚举a和b的匹配

2.打乱匹配顺序

假设我们已经完成了操作1，我们来计算每个匹配所能产生的期望合法方案

尝试转化一下模型：对于一个给定的匹配，我们从aiai向bibi连一条有向边，可以发现这个图最终由若干个环和若干条链构成，且链的顺序是唯一的

假设有e个Ai=Bi=1Ai=Bi=1，m个Ai=1,Bi=0Ai=1,Bi=0,可以发现边数为e+m，图由m条链和若干环组成

f[i][j]表示前i条链分到j个点的期望合法方案

那么有转移：f[i][j]=∑u≤ju=0f[i−1][j−u]/(u+1)
最终的答案即为e!∗m!∗(e+m)!∗∑j≤ej=0f[m][j]

代码

#include<bits/stdc++.h>
#define ll long long
using namespace std;
const int N=20077,p=998244353;
const ll top=1e17;
int power(int a,int b)
{
    int ass=1;
    for(;b;b>>=1,a=1ll*a*a%p) if(b&1) ass=1ll*ass*a%p;
    return ass;
}
char A[N],B[N];
int fac[N],unfac[N],val[N],l,n,m,c,ass;
int C(int n,int m)
{
    return 1ll*fac[n]*unfac[m]%p*unfac[n-m]%p;
}
int work(int n,int m,int k)
{
    return 1ll*C(m,k)*fac[k]%p*fac[m-k]%p*fac[m-k]%p*C(m+n,m-k)%p;
} 
int main()
{
    fac[0]=1;
    for(int i=1; i<=N; i++) fac[i]=1ll*fac[i-1]*i%p;
    unfac[N-1]=power(fac[N-1],p-2);
    for(int i=N-1;i;i--) unfac[i-1]=1ll*unfac[i]*i%p;
    scanf("%s%s",A+1,B+1);
    l=strlen(A+1);
    for(int i=1; i<=l; i++)
    {
        if(A[i]=='1'&&B[i]=='0') n++;
        if(A[i]=='1'&&B[i]=='1') m++;
    }
    for(int i=0; i<=n; i++) val[i]=1ll*C(n,i)*power(i,n)%p;
    for(int k=0; k<=m; k++)
    {
        ll tmp=0;
        for(int i=0; i<=n; i++)
        {
            if((n-i)&1) tmp-=val[i];else tmp+=val[i];
            if(tmp<-top||tmp>top) tmp%=p;
        }
        tmp=(tmp%p+p)%p;
        ass=(ass+1ll*fac[n]*tmp%p*work(n,m,k))%p;
        for(int i=0; i<=n; i++) val[i]=1ll*val[i]*i%p;
    }
    printf("%d\n",ass);
}