洛谷 P4705 玩游戏

洛谷 P4705 玩游戏


首先随便推一下式子:

\(ans_k=\sum_{i=1}^n\sum_{j=1}^m(a_i+b_j)^k\)
\(=\sum_{i=1}^n\sum_{j=1}^m\sum_{l=0}^ka_i^lb_j^{k-l}C_{k}^l\)
\(=\sum_{l=0}^kC_{k}^l(\sum_{i=1}^na_i^l)(\sum_{j=1}^mb_j^{k-l})\)

\(A_i=\sum_{j=1}^na_j^i,B_i=\sum_{j=1}^mb_j^i\)

\(ans_k=k!\sum_{l=0}^k\frac{A_{l}}{l!}\frac{B_{k-l}}{l!}\)

显然ntt就完事了,但是\(A,B\)怎么求。。然后膜题解去了

构造生成函数\(F(x)=A_0x^0+A_1x^1+A_2x^2+\cdots+A_{\infty}x^{\infty}\)

那么\(F(x)\)\(n\)个等比数列的和

\(F_i(x)=a_i^0x^0a_i^1x^1+a_i^2x^2+\cdots+a_i^{\infty}x^{\infty}\),显然\(F_i(x)=\frac 1{1-a_ix}\)

\(F(x)=\sum_{i=1}^nF_i(x)=\sum_{i=1}^n\frac 1{1-a_ix}\)

用分治NTT,每次合并两个,大力通分,最后分母乘分子逆元即是答案。。。注意项数要加到\(t\)

卡了一波常但常数依旧大。。

#include<bits/stdc++.h>
#define il inline
#define vd void
#define mod 998244353
#define G 3
#define iG 332748118
#define poly std::vector<int>
typedef long long ll;
il ll gi(){
    ll x=0,f=1;
    char ch=getchar();
    while(!isdigit(ch))f^=ch=='-',ch=getchar();
    while(isdigit(ch))x=x*10+ch-'0',ch=getchar();
    return f?x:-x;
}
il int pow(int x,int y){
    int ret=1;
    while(y){
        if(y&1)ret=1ll*ret*x%mod;
        x=1ll*x*x%mod;y>>=1;
    }
    return ret;
}
#define maxn 262147
int a[100010],b[100010];
poly pA,pB;
int rev[maxn],_lstN,P[maxn],iP[maxn];
il vd ntt(int*A,int N,int t){
    for(int i=0;i<N;++i)if(rev[i]>i)std::swap(A[i],A[rev[i]]);
    for(int o=1;o<N;o<<=1){
        int W=t?P[o]:iP[o];
        for(int*p=A;p!=A+N;p+=o<<1)
            for(int i=0,w=1;i<o;++i,w=1ll*w*W%mod){
                int t=1ll*w*p[i+o]%mod;
                p[i+o]=(p[i]-t+mod)%mod;p[i]=(p[i]+t)%mod;
            }
    }
    if(!t){
        int inv=pow(N,mod-2);
        for(int i=0;i<N;++i)A[i]=1ll*A[i]*inv%mod;
    }
}
int N,lg;
il vd setN(int n){
    N=1,lg=0;
    while(N<n)N<<=1,++lg;
    if(N!=_lstN)for(int i=0;i<N;++i)rev[i]=(rev[i>>1]>>1)|((i&1)<<lg-1);
}
il vd ntt(poly&a,int t){
    static int A[maxn];
    for(int i=0;i<a.size();++i)A[i]=a[i];memset(A+a.size(),0,4*(N-a.size()));
    ntt(A,N,t);
    a.resize(N);
    for(int i=0;i<N;++i)a[i]=A[i];
    int s=a.size();while(s&&!a[s-1])--s;
    a.resize(s);
}
il poly mul(poly a,poly b,int newn){
    setN(a.size()+b.size()-1);
    ntt(a,1),ntt(b,1);
    for(int i=0;i<N;++i)a[i]=1ll*a[i]*b[i]%mod;
    ntt(a,0);a.resize(newn);
    return a;
}
il poly operator+(poly a,const poly&b){
    if(a.size()<b.size())a.resize(b.size());
    for(int i=0;i<a.size();++i)if(i<b.size())a[i]=(a[i]+b[i])%mod;
    return a;
}
il poly operator-(poly a,const poly&b){
    if(a.size()<b.size())a.resize(b.size());
    for(int i=0;i<a.size();++i)if(i<b.size())a[i]=(a[i]-b[i]+mod)%mod;
    return a;
}
il poly operator*(poly a,int b){
    for(auto&i:a)i=1ll*i*b%mod;
    return a;
}
il poly qiudao(poly a){
    for(int i=0;i<a.size()-1;++i)a[i]=1ll*a[i+1]*(i+1)%mod;
    a.erase(a.end()-1);
    return a;
}
il poly jifen(poly a){
    a.insert(a.begin(),0);
    for(int i=1;i<a.size();++i)a[i]=1ll*a[i]*pow(i,mod-2)%mod;
    return a;
}
il poly getinv(poly a){
    if(a.size()==1)return poly(1,pow(a[0],mod-2));
    int n=a.size(),m=a.size()+1>>1;
    poly _a(m);
    for(int i=0;i<m;++i)_a[i]=a[i];
    poly b=getinv(_a);
    setN(n+m*2-2);
    ntt(a,1);ntt(b,1);
    for(int i=0;i<N;++i)a[i]=1ll*a[i]*b[i]%mod*b[i]%mod;
    ntt(a,0),ntt(b,0);
    a.resize(n);
    return b*2-a;
}
il poly getln(poly a){
    return jifen(mul(qiudao(a),getinv(a),a.size()));
}
il vd poly_init(){
    for(int i=1;i<maxn;i<<=1)P[i]=pow(G,(mod-1)/(i<<1)),iP[i]=pow(iG,(mod-1)/(i<<1));
}
il std::pair<poly,poly>divide(int l,int r,int*a){
    if(l==r){
        poly ret1,ret2;
        ret1.push_back(1);
        ret2.push_back(1);ret2.push_back(mod-a[l]);
        return {ret1,ret2};
    }
    int mid=(l+r)>>1;
    auto A=divide(l,mid,a),B=divide(mid+1,r,a);
    setN(A.second.size()+B.second.size()-1);
    ntt(A.first,1),ntt(B.first,1),ntt(A.second,1),ntt(B.second,1);
    for(int i=0;i<N;++i){
        int a=A.first[i],b=A.second[i],c=B.first[i],d=B.second[i];
        A.first[i]=(1ll*a*d+1ll*b*c)%mod;
        A.second[i]=1ll*b*d%mod;
    }
    ntt(A.first,0),ntt(A.second,0);
    return A;
}
int fact[100010],ifact[100010];
int main(){
#ifdef XZZSB
    freopen("in.in","r",stdin);
    freopen("out.out","w",stdout);
#endif
    poly_init();
    int n=gi(),m=gi();
    for(int i=1;i<=n;++i)a[i]=gi();
    for(int i=1;i<=m;++i)b[i]=gi();
    int t=gi()+1;
    auto _A=divide(1,n,a),_B=divide(1,m,b);
    _A.second.resize(t);_B.second.resize(t);
    poly A=mul(_A.first,getinv(_A.second),t),B=mul(_B.first,getinv(_B.second),t);
    fact[0]=1;for(int i=1;i<t;++i)fact[i]=1ll*fact[i-1]*i%mod;
    ifact[t-1]=pow(fact[t-1],mod-2);for(int i=t-2;~i;--i)ifact[i]=1ll*ifact[i+1]*(i+1)%mod;
    n=pow(n,mod-2),m=pow(m,mod-2);
    for(int i=0;i<A.size();++i)A[i]=1ll*A[i]*ifact[i]%mod*n%mod;
    for(int i=0;i<B.size();++i)B[i]=1ll*B[i]*ifact[i]%mod*m%mod;
    A=mul(A,B,t);
    for(int i=1;i<t;++i)printf("%d\n",1ll*A[i]*fact[i]%mod);
    return 0;
}
posted @ 2019-06-09 15:37  菜狗xzz  阅读(170)  评论(0编辑  收藏  举报