[多维FFT Bluestein′s Algorithm] Codechef October Challenge 2017 .Chef and Horcrux

题目里那个

其实不重要，只要能算出 pi 就行了

pi 的话发现就是个多维FFT的转移形式

到Hillan大佬博客里拷了个代码改一改就好了…

Bluestein算法里的FFT可以用暴力卷积代替，因为数据范围小FFT常数大

#include<cstdio>#include<iostream>#include<cstring>#include<cstdlib>#include<algorithm>#include<cmath>using namespace std;const int P=330301441;int rev[3200005];int w[2][3200005];typedef long long ll;inline int Pow(int x,ll y){    int ret=1;    for(;y;y>>=1,x=1LL*x*x%P) if(y&1) ret=1LL*ret*x%P;    return ret;}int c[3200005];int C,K,G,InvG;int tmp[3200005],tmp1[3200005],B1[3200005],B2[3200005],nxt[1110];int N,inv;inline void DFT(int* a,int d,int n,int fl){    int Nf=2*n,*A=tmp,*B=tmp1,*BB=B1;    if(fl==-1) swap(A,B),BB=B2;    for(int i=0,j=0;i<n;i++,j=nxt[i])             c[n-i]=1LL*a[i*d]*B[j]%P;    for(int i=n;i<2*n;i++){        int *pos=a+(i-n)*d;        *pos=0;        for(int j=0;j<=i;j++)            *pos=(*pos+1LL*BB[i-j]*c[j])%P;    }    for(int i=0,j=0;i<n;i++,j=nxt[i])             a[i*d]=1LL*a[i*d]*B[j]%P;    if(fl==-1)      for(int i=0;i<n;i++)a[i*d]=1LL*inv*a[i*d]%P;    for(int i=0;i<(1<<N);i++) c[i]=0;}void FWT(int* a,int n,int m,int fl){    int len=1;    for(int i=0;i<m;len*=K,i++)        for(int j=0;j<n;j+=len*K)            for(int k=0;k<len;k++)                DFT(a+j+k,len,K,fl);}int Da[3200005];int no[3200005];int num[2200010];int ans[3200005];inline char nc(){    static char buf[100000],*p1=buf,*p2=buf;    return p1==p2&&(p2=(p1=buf)+fread(buf,1,100000,stdin),p1==p2)?EOF:*p1++; }inline void rea(int &x){    char c=nc(); x=0;    for(;c>'9'||c<'0';c=nc());for(;c>='0'&&c<='9';x=x*10+c-'0',c=nc());}inline void rea(long long &x){    char c=nc(); x=0;    for(;c>'9'||c<'0';c=nc());for(;c>='0'&&c<='9';x=x*10+c-'0',c=nc());}int main(){    int T; rea(T);    while(T--)    {        int n,m=0; long long x;        rea(n); rea(K); rea(x);        for(int i=1,x;i<=n;i++)            rea(x),num[x]++;        ll pro=1,INV=Pow(Pow(2,n),P-2);        for(n=1;n<=100000;n*=K) m++;        for(int i=0;i<=n;i++) pro=1LL*pro*(num[i]=Pow(2,num[i]))%P;        for(int i=0,*j=Da;i<=n;i++,j++){            pro=pro*Pow(num[i],P-2)%P;            *j=pro*INV%P;            pro=pro*(num[i]-1)%P;            }            C=n*2;            G=Pow(22,(P-1)/2/K),InvG=Pow(G,P-2);        N=1; inv=Pow(K,P-2);        for(;(1<<N)<(3*K);N++);        int Nf=2*K,G=Pow(22,(P-1)/Nf),InvG=Pow(G,P-2);        tmp[0]=tmp1[0]=1;        for(int i=1;i<Nf;i++) tmp[i]=1LL*tmp[i-1]*G%P,tmp1[i]=1LL*tmp1[i-1]*InvG%P;            for(int i=0,j=0;i<2*K;i++,j=(j+(i<<1)+1)%Nf) B1[i]=tmp[j],B2[i]=tmp1[j];        for(int i=0,j=0;i<K;i++,j=nxt[i]) nxt[i+1]=((j+(i+1<<1)+1)%Nf);            FWT(Da,n,m,1);      for(int i=0;i<=n;i++) Da[i]=Pow(Da[i],x%(P-1));      FWT(Da,n,m,-1);      int ANS=0;      for(int i=0;i<=n;i++)        ANS=(ANS+1LL*Pow(i,2*i)*Pow(Da[i],3*i))%P;      cout<<ANS<<endl;      for(int i=0;i<=n;i++) num[i]=0;    }    return 0;}