2020 ccpc 网络赛 1013 Residual Polynomial
题意:
为什么人均都会 NTT+分治啊,考试最后半个小时推出来式子,结果打不完了。
通过仔细分析或者手动模拟,我们可以发现,a[ i ]对于 f[ n ] [ j ] (i>=j) 的系数的贡献是由 b[p1]*b[p2]*b[p3]……*b[p(i-j)]*c[q1]*c[q2]*……*c[q(n-1-i-j)]*a[i]*( i*(i-1)*……*(j+1) ),其中q与p中的元素各不相同,可以从类似组合数的思想理解它的贡献方式。
不妨设H[x]为 i-j==x 时a[i]前的那一大串系数。我们会发现,f[ n ][ i ]* i! =H[ j ]*a[ i+j ]*(i+j)! (0<=j<=n-i),我们将H[j]与H[n-j]调换,式子就变成了 f[n][i]*i!=H[n-j]*a[i+j]*(i+j)!,我们会发现这个式子可以被NTT卷起来求值,我们剩下的问题就是如何求H。
暴力求H的方法是用背包,复杂度n^2,显然不太优美,因此我们可以用分治+NTT的方法,将复杂度降下去,小范围直接暴力跑就可以了。
1 #include <bits/stdc++.h> 2 #define N 800005 3 using namespace std; 4 inline int read(){ 5 int x=0,f=1;char ch=getchar(); 6 while(ch<'0'||ch>'9'){if(ch=='-')f=-1;ch=getchar();} 7 while(ch>='0'&&ch<='9'){x=x*10+ch-'0';ch=getchar();} 8 return x*f; 9 } 10 int T,n; 11 int a[N],b[N],c[N]; 12 int rev[N]; 13 long long A[N],B[N]; 14 int F[2][18][N/4]; 15 const int p=998244353; 16 long long ksm(long long x,long long z) 17 { 18 if(!x)return 1; 19 long long ans=1; 20 while(z>0) 21 { 22 if(z&1) 23 { 24 ans*=x; 25 ans%=p; 26 } 27 x*=x; 28 x%=p; 29 z>>=1; 30 } 31 return ans; 32 } 33 void NTT(long long *A,int n,int op) 34 { 35 for(int i=1;i<n;i++) 36 { 37 if(rev[i]>i) swap(A[i],A[rev[i]]); 38 } 39 for(int k=2;k<=n;k<<=1) 40 { 41 long long wn=ksm(3,op==1?(p-1)/k:(p-1)-(p-1)/k); 42 for(int j=0;j<n;j+=k) 43 { 44 long long w=1,x,y; 45 for(int i=0;i<(k>>1);i++) 46 { 47 x=A[i+j],y=w*A[i+j+(k>>1)]%p; 48 A[i+j]=(x+y)%p; 49 A[i+j+(k>>1)]=(x-y+p)%p; 50 w*=wn,w%=p; 51 } 52 } 53 } 54 if(op==-1) 55 { 56 int inv=ksm(n,p-2); 57 for(int i=0;i<n;i++) A[i]=1ll*A[i]*inv%p; 58 } 59 } 60 void work(int *AA,int *BB,int n) 61 { 62 int t,l; 63 for(t=1,l=0;t<=n;t<<=1,l++); 64 l--; 65 for(int i=1;i<t;i++)rev[i]=(rev[i>>1]>>1)|((1&i)<<l); 66 for(int i=0;i<=n;i++)A[i]=AA[i],B[i]=BB[i]; 67 for(int i=n+1;i<=t;i++) A[i]=B[i]=0; 68 NTT(A,t,1); 69 NTT(B,t,1); 70 for(int i=0;i<t;i++) A[i]=A[i]*B[i]%p; 71 NTT(A,t,-1); 72 } 73 void solve(int l,int r,int op,int dep) 74 { 75 if(r-l+1<=100) 76 { 77 F[op][dep][0]=1; 78 for(int i=l;i<=r;i++) 79 { 80 F[op][dep][i-l+1]=0; 81 for(int j=i-l+1;j;j--) 82 { 83 F[op][dep][j]=(1ll*F[op][dep][j-1]*b[i]%p+1ll*F[op][dep][j]*c[i]%p)%p; 84 } 85 F[op][dep][0]=1ll*F[op][dep][0]*c[i]%p; 86 } 87 F[op][dep][r-l+2]=F[op][dep][r-l+3]=0; 88 return; 89 } 90 int mid=(l+r)>>1; 91 solve(l,mid,0,dep+1); 92 solve(mid+1,r,1,dep+1); 93 work(F[0][dep+1],F[1][dep+1],r-l+1); 94 95 F[0][dep+1][0]=F[1][dep+1][0]=0; 96 for(int i=0;i<=r-l+1;i++) F[op][dep][i]=A[i],A[i]=B[i]=F[0][dep+1][i]=F[1][dep+1][i]=0; 97 A[0]=0; 98 99 F[op][dep][r-l+2]=F[op][dep][r-l+3]=0; 100 } 101 int G[N],jc[N],C[N]; 102 int main(){ 103 // freopen("test.in","r",stdin); 104 // freopen("1.out","w",stdout); 105 // scanf("%d",&T); 106 T=read(); 107 jc[0]=1; 108 for(int i=1;i<=100000;i++) jc[i]=1ll*jc[i-1]*i%p; 109 while(T--) 110 { 111 //scanf("%d",&n); 112 n=read(); 113 for(int i=0;i<=n;i++) 114 { 115 //scanf("%d",&a[i]); 116 a[i]=read(); 117 } 118 for(int i=1;i<=n-1;i++) 119 { 120 //scanf("%d",&b[i]); 121 b[i]=read(); 122 } 123 for(int i=1;i<n;i++) 124 { 125 //scanf("%d",&c[i]); 126 c[i]=read(); 127 } 128 memset(G,0,sizeof(G)); 129 memset(F,0,sizeof(F)); 130 memset(C,0,sizeof(C)); 131 memset(A,0,sizeof(A)); 132 memset(B,0,sizeof(B)); 133 long long sm=1;; 134 solve(1,n-1,0,0); 135 for(int i=0;i<n;i++) G[i]=F[0][0][i]; 136 // for(int i=0;i<n;i++) cout<<G[i]<<' '; 137 // cout<<endl; 138 for(int i=0;i<=n/2;i++) swap(G[i],G[n-i]); 139 for(int i=0;i<=n;i++) C[i]=1ll*jc[i]*a[i]%p; 140 work(G,C,n*2); 141 for(int i=n;i<=n*2;i++) 142 { 143 printf("%lld",A[i]*ksm(jc[i-n],p-2)%p); 144 if(i!=n*2) printf(" "); 145 } 146 printf("\n"); 147 148 } 149 return 0; 150 }