hdu 6036 NTT

来源：互联网发布：samba windows共享编辑：程序博客网时间：2024/06/05 03:18
#pragma comment(linker, "/STACK:102400000,102400000")#include<bits/stdc++.h>using namespace std;typedef long long ll;const int maxn=2*1e5+5;const ll mod=985661441;const ll P=mod;//费马素数const ll G=3;//原根ll inv[maxn],facinv[maxn],fac[maxn],e[maxn],g[maxn],ans[15],pa[2][maxn];ll a[2*maxn],b[2*maxn],wn[26];ll mul(ll x,ll y){    return (x%P)*(y%P)%P;}ll quick_mod(ll a,ll b){  ll res=1;  while(b)  {    if(b&1)        res=res*a%P;    a=a*a%mod;    b>>=1;  }  return res;}void GetWn(){    for(int i=0;i<26;i++)    {        int t =1<<i;        wn[i]=quick_mod(G,(P-1)/t);    }}void Rader(ll a[],int len){    int j=len>>1;    for(int i=1;i<len-1;i++)    {        if(i<j) swap(a[i],a[j]);        int k=len>>1;        while(j>=k)        {            j-=k;            k>>=1;        }        if(j<k)            j+=k;    }}void NTT(ll a[],int len,int on){    Rader(a,len);    int id=0;    for(int h=2;h<=len;h<<=1)    {        id++;        for(int j=0;j<len;j+=h)        {            ll w=1;            for(int k=j;k<j+h/2;k++)            {                ll u=a[k];                ll t=mul(w,a[k+h/2]);                a[k]=(u+t)%P;                a[k+h/2]=((u-t)%P+P)%P;                w=mul(w,wn[id]);            }        }    }    if(on==-1)    {        for(int i=1;i<len/2;i++)            swap(a[i],a[len-i]);        ll Inv=quick_mod(len,P-2);        for(int i=0;i<len;i++)            a[i]=mul(a[i],Inv);    }}void Conv(ll a[],ll b[],int n){    NTT(a,n,1);    NTT(b,n,1);    for(int i=0;i<n;i++)        a[i]=mul(a[i],b[i]);    NTT(a,n,-1);}void init(){    fac[0]=facinv[0]=fac[1]=facinv[1]=inv[1]=inv[0]=1;    for(int i=2;i<maxn;i++)    {        inv[i]=(mod-mod/i)*inv[mod%i]%mod;        fac[i]=i*fac[i-1]%mod;        facinv[i]=inv[i]*facinv[i-1]%mod;    }}ll cob(ll n,ll m){    if(n<0||m<0)        return 0;    return fac[n]*facinv[n-m]%mod*facinv[m]%mod;}int main(){    ///////////////////////////////////////////////////    GetWn();    init();    //////////////////////////////////////////////////    int m,k,i,j,cas=0;    ll n;    while(scanf("%d%d",&m,&k)!=EOF)    {        for(i=n=0;i<m;i++)        {            scanf("%*lld%lld",e+i);            n+=e[i];        }        for(i=0;i<=n;i++)        {            g[i]=1;            for(j=0;j<m;j++)                g[i]=g[i]*cob(e[j]+i-1,(ll)i-1)%mod;        }        //////////////////////////////////////////////////        int len=1;        while(len<2*(n+1))            len<<=1;        //////////////////////////////////////////////////        for(i=0;i<=n;i++)        {            a[i]=(i%2?mod-facinv[i]:facinv[i]);            b[i]=g[i]*facinv[i]%mod;        }        //////////////////////////////////////////////////        for(i=n+1;i<len;i++)            a[i]=b[i]=0;        Conv(a,b,len);        //////////////////////////////////////////////////存在a数组中        for(i=0;i<=n;i++)            a[i]=a[i]*fac[i]%mod;        memset(ans,0,sizeof(ans));        a[n+1]=0;        int pre=1,cur=0;        for(i=1;i<=n;i++)        {            pre^=1,cur^=1;            pa[cur][0]=1;            for(j=1;j<=k;++j)                pa[cur][j]=pa[cur][j-1]*a[i]%mod;            if(i>1)                for(j=1;j<=k;++j)                    ans[j]=(ans[j]+pa[cur][j-1]*pa[pre][k-j+1]%mod)%mod;        }        ans[1]=(ans[1]+pa[cur][k])%mod;        printf("Case #%d: ",++cas);        for(i=1;i<k;i++)            printf("%d ",ans[i]);        printf("%d\n", ans[i]);    }    return 0;}
阅读全文
0 0