题目里那个
其实不重要,只要能算出 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;
}