51nod1172 Partial Sums V2

推一下式子发现是裸的FFT，$ans[k]=sum_{i}sum_{j}[i+j=k]a[i]*C_{m-1+j}^{j}$

比较坑爹的就是这个模数，于是我们上任意模数的FFT

任意模数的FFT目的就是降低卷积中的元素上界，我们设$P=lfloor sqrt{mod} floor$，

我们将原函数中的系数变成两个$a1[i]=a[i]/P,a2[i]=a[i]%P,b1[i]=b[i]/P,b2[i]=b[i]%P$

这样我们新卷积出来的值的上限是$n*mod$，

$P^2$的系数是$a1*b1$，$P$的系数是$a1*b2+a2*b1$，$1$的系数是$a2*b2$，然后我们再分别卷积，最后统计答案即可

一开始炸精了。。改成预处理单位复数根就好了。

#include <cstdio>
#include <cstring>
#include <iostream>
#include <algorithm>
#include <cmath>
#define mod 1000000007
#define MAXN 131072
using namespace std;
const double pi=acos(-1.0);
const int P=31622;
struct cmx{
    double x,y;
    cmx(){}
    cmx(double a,double b){x=a,y=b;}
    cmx operator + (cmx a){return cmx(x+a.x,y+a.y);}
    cmx operator - (cmx a){return cmx(x-a.x,y-a.y);}
    cmx operator * (cmx a){return cmx(x*a.x-y*a.y,x*a.y+y*a.x);}
    cmx operator / (double a){return cmx(x/a,y/a);}
}A[MAXN],B[MAXN],C[MAXN],D[MAXN],E[MAXN],F[MAXN],G[MAXN],wn[MAXN],inv[MAXN];
int n,m,N,rev[MAXN];
long long ans[MAXN],a[MAXN],b[MAXN];
void fft(cmx *a,int o, cmx *wn){
    register int i,j,k;
    cmx t;
    for(i=0;i<N;i++)if(i<rev[i])swap(a[i],a[rev[i]]);
    for(k=2;k<=N;k<<=1)
        for(j=0;j<N;j+=k)
            for(i=0;i<(k>>1);i++){
                t=wn[N/k*i]*a[i+j+(k>>1)];
                a[i+j+(k>>1)]=a[i+j]-t;
                a[i+j]=a[i+j]+t;
            }
    if(o==-1)for(i=0;i<N;i++)a[i]=a[i]/N;
}
void FFT(long long *a,long long *b,long long *ans){
    register int i;
    for(i=0;i<N;i++){
        A[i]=cmx(a[i]/P,0);
        B[i]=cmx(a[i]%P,0);
        C[i]=cmx(b[i]/P,0);
        D[i]=cmx(b[i]%P,0);
    }
    fft(A,1,wn);fft(B,1,wn);fft(C,1,wn);fft(D,1,wn);
    for(i=0;i<N;i++){
        E[i]=A[i]*C[i];
        F[i]=A[i]*D[i]+B[i]*C[i];
        G[i]=B[i]*D[i];
    }
    fft(E,-1,inv);fft(F,-1,inv);fft(G,-1,inv);
    register long long tmp;
    for(i=0;i<n;i++){
        tmp=1ll*round(E[i].x);
        ans[i]=tmp%mod*P%mod*P%mod;
        tmp=1ll*round(F[i].x);
        (ans[i]+=tmp%mod*P%mod)%=mod;
        (ans[i]+=1ll*round(G[i].x))%=mod;
    }
}
long long qp(long long a,int b){
    long long c=1;
    while(b){
        if(b&1)c=c*a%mod;
        a=a*a%mod; b>>=1;
    }return c;
}
int main(){
    scanf("%d%d",&n,&m);
    register int i;
    for(i=0;i<n;i++)scanf("%lld",&a[i]);
    b[0]=1;
    for(i=1;i<n;i++)
        b[i]=1ll*b[i-1]*(m-1+i)%mod*qp(i,mod-2)%mod;
    for(N=1;N<(n<<1);N<<=1);
    for(i=0;i<N;i++){
        if(i&1)rev[i]=(rev[i>>1]>>1)|(N>>1);
        else rev[i]=rev[i>>1]>>1;
    }
    for(i=0;i<N;i++)
        wn[i]=cmx(cos(2*pi/N*i),sin(2*pi/N*i)),
        inv[i]=cmx(cos(2*pi/N*i),-sin(2*pi/N*i));
    FFT(a,b,ans);
    for(i=0;i<n;i++)printf("%lld
",ans[i]);
    return 0;
}