2019雅礼集训 D5T3 permutation [DP,数学,FFT]

题目描述:

2019雅礼集训 D5T3 permutation [DP,数学,FFT]

样例:

input1:
1 1 1

output1:
1

input2:
2 1 1

output2:
0

input3:
2 2 1

output3:
1

input4:
5 2 2

output4:
22

数据范围与约定:

2019雅礼集训 D5T3 permutation [DP,数学,FFT]


又是原题……CF960G

看到最后一个子任务了吗?毒瘤出题人卡常数……

此题首先注意到有一个子任务是\(a+b=n+1\),思考什么时候才能出现这种情况?

最大的放在中间,左边a-1个单调递增,右边b-1个单调递减。

那么这个子任务的答案就出来了:\(C(n-1,a-1)\)

同时,它也引导我们研究最大的元素的摆放。

考虑DP:设\(dp_{i,j}\)表示从左到右有i个jm数,共用了j个数,此时的方案数。(不要问我为什么状态似乎i,j反了,这是为了接下来好看作准备)

那么\(ans=\sum_{i=0}^{n-1} dp_{a-1,i}\cdot dp_{n-i-1,b-1}\cdot{{n-1}\choose i}\),即枚举最大数的位置。

考虑向原数列插入一个最大的元素\(a_{max}\):枚举它左边有多少个元素,那么\(a_{max}\)必定会对jm数的个数作出1的贡献,右边没有,全靠左边。

于是有了递推公式:\(dp_{i,j}=\sum_{k=0}^{j-1}({{j-1}\choose k}\cdot dp_{i-1,k}\cdot (j-k-1)!)\)

可以看出这是一个优秀的\(O(n^3)\)的DP方程,显然会挂。

优化一波:把组合数拆开,可以得到\(dp_{i,j}=(j-1)!\sum_{k=0}^{j-1}\frac{ dp_{i-1,k}}{k!}\)

此时记个前缀和,就可以\(n^2\)递推了。结合数据分治,能得到不错的分数。

然而这个优化对得到正解一点用也没有

优化之路受阻,我们尝试将dp数组打印出来找规律。

2019雅礼集训 D5T3 permutation [DP,数学,FFT]

如果你经验丰富,就可以看出:这不正是第一类斯特林数吗?!可以\(O(n\log n)\)快速求出一整行的啊!!

然而,答案的统计需要一整列……

经验丰富的选手可以再次看出:答案的统计是可以优化的!

2019雅礼集训 D5T3 permutation [DP,数学,FFT]

2019雅礼集训 D5T3 permutation [DP,数学,FFT]

于是我们只需要求出一个斯特林数就可以啦!

做完啦!!完结撒花!!

并没有。

对于经验丰富的选手,可能的确是做完了。但像我这种菜鸡仍然在懵逼:蛤?咋求啊?

网上一搜:第一类斯特林数:

一、分治FFT瞎搞,\(O(n\log ^2 n)\)

然而我并没有看懂他们是怎么瞎搞的……知道的人麻烦在评论区回复一下,谢谢。

二、倍增+FFT,\(O(n\log n)\)

终于看懂了,但是好难写啊……

然而有些晚了,我不想写这式子怎么推了……自行百度吧……

代码:

#include<bits/stdc++.h>
namespace my_std{
    using namespace std;
    #define mod 998244353
    #define pii pair<int,int>
    #define fir first
    #define sec second
    #define MP make_pair
    #define rep(i,x,y) for (int i=(x);i<=(y);i++)
    #define drep(i,x,y) for (int i=(x);i>=(y);i--)
    #define go(x) for (int i=head[x];i;i=edge[i].nxt)
    #define sz 500500
    typedef long long ll;
    template<typename T>
    inline void read(T& t)
    {
        t=0;char f=0,ch=getchar();
        double d=0.1;
        while(ch>'9'||ch<'0') f|=(ch=='-'),ch=getchar();
        while(ch<='9'&&ch>='0') t=t*10+ch-48,ch=getchar();
        if(ch=='.')
        {
            ch=getchar();
            while(ch<='9'&&ch>='0') t+=d*(ch^48),d*=0.1,ch=getchar();
        }
        t=(f?-t:t);
    }
    template<typename T,typename... Args>
    inline void read(T& t,Args&... args){read(t); read(args...);}
    void file()
    {
        #ifndef ONLINE_JUDGE
        freopen("a.txt","r",stdin);
        #endif
    }
    inline ll mul(ll a,ll b){ll d=(ll)(a*(double)b/mod+0.5);ll ret=a*b-d*mod;if (ret<0) ret+=mod;return ret;}
}
using namespace my_std;

#define int ll

ll fac[sz],_fac[sz];

ll ksm(ll x,int y)
{
    ll ret;
    for (ret=1;y;y>>=1,x=x*x%mod) if (y&1) ret=ret*x%mod;
    return ret;
}
ll inv(ll x){return ksm(x,mod-2);}

void init()
{
    fac[0]=_fac[0]=1;
    rep(i,1,sz-5) _fac[i]=inv(fac[i]=fac[i-1]*i%mod);
}

const ll g=3;
int r[sz],limit;
void NTT(ll *a,int type)
{
    rep(i,0,limit-1) if (r[i]<i) swap(a[i],a[r[i]]);
    for (int mid=1;mid<limit;mid<<=1)
    {
        ll Wn=ksm(g,(mod-1)/mid>>1);if (type==-1) Wn=inv(Wn);
        for (int len=mid<<1,j=0;j<limit;j+=len)
        {
            ll w=1;
            for (int k=0;k<mid;k++,w=1ll*w*Wn%mod)
            {
                ll x=a[j+k],y=1ll*w*a[mid+j+k];
                a[j+k]=(x+y+mod)%mod;a[j+k+mid]=(x-y+mod)%mod;
            }
        }
    }
    if (type==-1)
    {
        ll t=inv(limit);
        rep(i,0,limit-1) a[i]=(1ll*a[i]*t%mod+mod)%mod;
    }
}
void NTT(ll *a,ll *b,int len)
{
    limit=1;int l=-1;
    while (limit<=len*2) limit<<=1,++l;
    rep(i,1,limit-1) r[i]=(r[i>>1]>>1)|((i&1)<<l);
    rep(i,len,limit) a[i]=b[i]=0;
    NTT(a,1);NTT(b,1);
    rep(i,0,limit) a[i]=1ll*a[i]*b[i]%mod;
    NTT(a,-1);
}

ll X[sz],a[sz],b[sz];
void work(int len)
{
    if (len==1) return (void)(X[1]=1);
    if (len&1)
    {
        work(len-1);
        X[len]=0;
        drep(i,len,1) (X[i]=X[i-1]+1ll*X[i]*(len-1)%mod)%=mod;
        return;
    }
    int n=len>>1;
    work(n);
    for (ll i=n,w=1;i>=0;i--,w=1ll*w*n%mod) a[i]=1ll*w*_fac[n-i]%mod;
    rep(i,0,n) b[i]=fac[i]*X[i]%mod;
    NTT(a,b,n+1);
    rep(i,0,n) a[i]=1ll*a[i+n]*_fac[i]%mod;
    NTT(X,a,n+1);
}

int n,A,B;

signed main()
{
    file();
    init();
    read(n,A,B);
    if (n==1) return puts(A==B&&A==1?"1":"0"),0;
    work(n-1);
    cout<<1ll*X[A+B-2]*fac[A+B-2]%mod*_fac[A-1]%mod*_fac[B-1]%mod;
}
上一篇:【jxoi2018】游戏 组合数学


下一篇:ERROR: ORA-28009: connection as SYS should be as SYSDBA or SYSOPER