zoukankan      html  css  js  c++  java
  • [题解] LuoguP5396 第二类斯特林数·列

    https://www.luogu.com.cn/problem/P5396

    Description:

    让你求第二类斯特林数第(k)列的前(n)

    (167772161)取模

    Solution

    考虑第(k)列斯特林数的生成函数

    [F_k(x) = sumlimits_{i} egin{Bmatrix} i \ k end{Bmatrix} x^i ]

    众所周知

    [egin{Bmatrix} i \ kend{Bmatrix} = kegin{Bmatrix} i-1 \ kend{Bmatrix} + egin{Bmatrix} i -1\ k-1end{Bmatrix} ]

    于是

    [egin{aligned}F_k(x) &= sumlimits_{i} left( kegin{Bmatrix} i-1\kend{Bmatrix} + egin{Bmatrix} i-1 \ k-1end{Bmatrix} ight)x^i \ &= left(sumlimits_{i} kegin{Bmatrix} i-1 \ kend{Bmatrix}x^i ight) + left( sumlimits_{i} egin{Bmatrix} i-1 \ k-1end{Bmatrix}x^i ight)\&= left(kxsumlimits_{i} egin{Bmatrix} i \ kend{Bmatrix}x^i ight) + left(xsumlimits_{i}egin{Bmatrix} i\ k-1end{Bmatrix}x^i ight)\&= kxF_k(x) + xF_{k-1}(x)end{aligned} ]

    移项后可以得到

    [F_k(x) = frac{xF_{k-1}(x)}{1-kx} ]

    边界就是(F_0(x) = 1)

    于是把上面那玩意儿一直展开下去,可以得到

    [F_k(x) = frac{x^k}{prod_{i=1}^{k} (1-ix)} ]

    分治乘搞出分母,然后求逆在乘上(x^k)就好了

    终于不卡分治乘了(泪目

    (Code)

    #include <bits/stdc++.h>
    using namespace std;
    const int N=1111111,MOD=167772161;
    inline int qpow(int a,int b=MOD-2,int m=MOD)
    {
        int ans=1%m;
        for(;b;b>>=1,a=1ll*a*a%m)
            if(b&1) ans=1ll*ans*a%m;
        return ans;
    }
    const int gn=3,ign=qpow(gn);
    namespace FFT {int r[N];}
    inline int glim(int n)
    {
        int lim=1;
        while(lim<=n) lim<<=1;
        return lim;
    }
    void fftinit(int n)
    {
        using namespace FFT;
        for(int i=0;i<n;i++) r[i]=r[i>>1]>>1|((i&1)?n>>1:0);
    }
    void fft(int *f,int n,int flg)
    {
        using namespace FFT;
        for(int i=0;i<n;i++) if(r[i]<i) swap(f[i],f[r[i]]);
        for(int len=2,k=1;len<=n;len<<=1,k<<=1)
        {
            int wn=qpow(flg==1?gn:ign,(MOD-1)/len);
            for(int i=0;i<n;i+=len)
                for(int j=i,w=1;j<i+k;j++,w=1ll*w*wn%MOD)
                {
                    int tmp=1ll*f[j+k]*w%MOD;
                    f[j+k]=f[j]-tmp<0?f[j]-tmp+MOD:f[j]-tmp;
                    f[j]=f[j]+tmp>=MOD?f[j]+tmp-MOD:f[j]+tmp;
                }
        }
        if(flg==-1)
        {
            int inv=qpow(n);
            for(int i=0;i<n;i++) f[i]=1ll*f[i]*inv%MOD;
        }
    }
    void ginv(const int *f,int n,int *g)
    {
        if(n==1){g[0]=qpow(f[0]);return;}
        ginv(f,(n+1)>>1,g);
        static int f1[N]; int lim=glim((n-1)<<1);
        fftinit(lim);
        for(int i=0;i<lim;i++) g[i]=i<n?g[i]:0,f1[i]=i<n?f[i]:0;
        fft(g,lim,1),fft(f1,lim,1);
        for(int i=0;i<lim;i++) g[i]=1ll*g[i]*(2ll-1ll*f1[i]*g[i]%MOD+MOD)%MOD;
        fft(g,lim,-1);
        for(int i=n;i<lim;i++) g[i]=0;
    }
    namespace PolyMul {int f[N],g[N];}
    struct poly
    {
        vector<int> a;
        inline int size() const {return a.size();}
        inline int deg() const {return size()-1;}
        inline int operator [] (int i) const {return i<size()?a[i]:0;}
        inline int &operator [] (int i) {assert(i<size());return a[i];}
        poly(int n=0,int v=0) {a=vector<int>(n,v);}
        void resize(int n) {a.resize(n);}
        void debug() {printf("Poly Debug: ");for(int i=0;i<size();i++) printf("%d%c",a[i]," 
    "[i==size()-1]);}
        void reverse() {for(int i=0,j=size()-1;i<j;i++,j--) swap(a[i],a[j]);}
        poly inv(int n)
        {
            static int f[N],g[N];
            memset(f,0,sizeof(f)),memset(g,0,sizeof(g));
            for(int i=0;i<n;i++) f[i]=i<size()?a[i]:0; // !!!!!! [] !!!!!!
            ginv(f,n,g); poly c(n);
            for(int i=0;i<n;i++) c[i]=g[i];
            return c;
        }
    };
    poly operator * (const poly &a,const poly &b)
    {
        using namespace PolyMul;
        int n=a.deg(),m=b.deg(),lim=glim(n+m);
        fftinit(lim);
        for(int i=0;i<lim;i++) f[i]=a[i],g[i]=b[i];
        fft(f,lim,1),fft(g,lim,1);
        for(int i=0;i<lim;i++) f[i]=1ll*f[i]*g[i]%MOD;
        fft(f,lim,-1);
        poly c(n+m+1);
        for(int i=0;i<c.size();i++) c[i]=f[i];
        return c;
    }
    poly solve(int l,int r)
    {
        if(l==r)
        {
            poly ans(2);
            ans[0]=1,ans[1]=MOD-l;
            return ans;
        }
        int mid=(l+r)>>1;
        return solve(l,mid)*solve(mid+1,r);
    }
    poly ans;
    int main()
    {
    #ifdef LOCAL
        freopen("a.in","r",stdin);
        freopen("a.out","w",stdout);
    #endif
        int n,k;scanf("%d%d",&n,&k);
        if(n-k+1<=0)
        {
            for(int i=0;i<=n;i++) printf("0 ");
            return 0;
        }
        ans=solve(1,k).inv(n-k+1);
        for(int i=0;i<k;i++) printf("0 ");
        for(int i=k;i<=n;i++) printf("%d ",ans[i-k]);
        return 0;
    }
    
  • 相关阅读:
    【tarjan】【树的直径】【CF】K. Königsberg Bridges
    【组合数学】【恒等式】简单和、交错和
    【组合数学】【恒等式】$sum_{k=0}^{r}C_m^k imes C_{n}^{r-k}=C_{m+n}^r$
    【组合数学】【恒等式】$C_{n}^{r} imes C_{n-r}^{k-r}=C_{n}^{k} imes C_k^{r}$
    【树形DP】D. Serval and Rooted Tree
    【图论】图的欧拉定理
    【图论】网络流解决二分图最大匹配量问题
    【计算几何】atan2函数
    【单峰计数DP】Problem F – Fabricating Sculptures
    Java基础之:自定义泛型
  • 原文地址:https://www.cnblogs.com/wxq1229/p/12952051.html
Copyright © 2011-2022 走看看