多项式
不知不觉,就开始学习多项式了,好快啊
重新拾起我之前学的百嘛不是的\(FFT\),新学\(NTT/FWT\)
开始吧学习博客
FFT:快速傅里叶变换
总的来说就是快速求取两个多项式的乘积
把板子粘在这里Luogu3803 多项式乘法
code
#include<bits/stdc++.h>
using namespace std;
#define fo(i,x,y) for(int i=(x);i<=(y);i++)
#define fu(i,x,y) for(int i=(x);i>=(y);i--)
const double pi=acos(-1.0);
const int N=1<<21;
struct coex{
double r,i;
coex(){}
coex(double x,double y){r=x;i=y;}
coex operator + (coex a){return coex(r+a.r,i+a.i);}
coex operator - (coex a){return coex(r-a.r,i-a.i);}
coex operator * (coex a){return coex(r*a.r-i*a.i,r*a.i+i*a.r);}
}a[N],b[N],w[N];
int af[N],la,lb,lim,len;
void fft(coex *a,int lim){
fo(i,0,lim-1)if(af[i]>i)swap(a[i],a[af[i]]);
for(int t=lim>>1,d=1;d<lim;d<<=1,t>>=1)
for(int i=0;i<lim;i+=(d<<1))
fo(j,0,d-1){
coex tmp=w[t*j]*a[i+j+d];
a[i+j+d]=a[i+j]-tmp;
a[i+j]=a[i+j]+tmp;
}
}
signed main(){
scanf("%d%d",&la,&lb);
fo(i,0,la)scanf("%lf",&a[i].r);
fo(i,0,lb)scanf("%lf",&b[i].r);
for(lim=1,len=0;lim<=la+lb;lim<<=1,len++);
fo(i,0,lim-1){
af[i]=(af[i>>1]>>1)|((i&1)<<(len-1));
w[i]=coex(cos(2.0*i*pi/lim),sin(2.0*i*pi/lim));
}
fft(a,lim);fft(b,lim);
fo(i,0,lim-1)a[i]=a[i]*b[i],w[i].i=-w[i].i;
fft(a,lim);
fo(i,0,la+lb)printf("%d ",(int)(a[i].r/lim+0.5));
}
需要注意这么几点:
1、\(lim\)要大于多项式次数,严格大于
2、搞清楚\(d\)和\(t\)的关系,以及各种2倍
3、背板子,背就完事了。
4、最后答案别忘了除以\(lim\),预处理单位根别忘了\(2*i*pi/lim\),别忘了除以\(lim\)
NTT:快速数论变换
code
void ntt(int *a,int lim){
fo(i,0,lim-1)if(af[i]>i)swap(a[i],a[af[i]]);
for(int t=lim>>1,d=1;d<lim;d<<=1,t>>=1)
for(int i=0;i<lim;i+=(d<<1))
fo(j,0,d-1){
int tmp=g[t*j]*a[i+j+d]%mod;
a[i+j+d]=(a[i+j]-tmp+mod)%mod;
a[i+j]=(a[i+j]+tmp)%mod;
}
}
int bas[M],ans[M];
int a[M],b[M];
void mul(int *x,int *y,int *z){
fo(i,0,lim-1)a[i]=y[i];
fo(i,0,lim-1)b[i]=z[i];
g[0]=1;g[1]=ksm(3,(mod-1)/lim,mod);
fo(i,2,lim-1)g[i]=g[i-1]*g[1]%mod;
ntt(a,lim);ntt(b,lim);
g[0]=1;g[1]=ksm(g[1],mod-2,mod);
fo(i,2,lim-1)g[i]=g[i-1]*g[1]%mod;
fo(i,0,lim-1)a[i]=a[i]*b[i]%mod;
ntt(a,lim);int inv=ksm(lim,mod-2,mod);
fo(i,0,lim-1)a[i]=a[i]*inv%mod;
fo(i,0,m-2)a[i]=(a[i]+a[i+m-1])%mod;
fo(i,0,m-2)x[i]=a[i];
}
注意:
1、注意模数,以及自己找原根
2、注意最后出来的时候除以\(lim\)
3、别忘了变换原根
FWT:快速沃尔什变换
发现我们所做的卷积都是\(c=\sum_{i=1}^{n}\sum_{j=1}^{i}a_j*b_{i-j}\)
也就是\(FFT\)和\(NTT\)能处理的加法卷积,如果是与?或?异或呢?
就要用到\(FWT\),也就是快速位运算卷积
|:或卷积
这个要是直接卷的话,好像有点困难,但是高维前缀和就是或卷积,所以我们可以用高维前缀和解决或卷积的问题
然而高维前缀和还是没有\(FWT\)来的直接
\(FWT\)也和\(FFT\)一样有变换和逆变换,这个直接看最上面我给出的那个博客吧
这次变换不能叫点值表达式了,应该叫它子集和表达式
相当于是枚举每一个二进制位来求子集和,逆变换拆回去的时候是一样的
code
void ort(int *a,int lim,int tp){//tp=1表示顺变换,tp=-1表示逆变换
for(int d=1;d<lim;d<<=1)
for(int i=0;i<lim;i+=(d<<1))
fo(j,0,d-1){
a[i+j+d]=(a[i+j+d]+a[i+j]*tp)%mod;
}
}
&:与卷积
这应该叫做高维后缀和吧,我们可以用高维后缀和解决与卷积问题
变换就是后缀和,逆变换一样
code
void andt(int *a,int lim,int tp){//tp=1表示顺变换,tp=-1表示逆变换
for(int d=1;d<lim;d<<=1)
for(int i=0;i<lim;i+=(d<<1))
fo(j,0,d-1){
a[i+j]=(a[i+j]+a[i+j+d]*tp)%mod;
}
}
^:异或卷积
这个应该叫高维后缀差分
好像有点难以概括,额,不概括了,直接看码吧
code
void xort(int *a,int lim,int tp){//tp=1表示顺变换,tp=0.5表示逆变换
for(int d=1;d<lim;d<<=1)
for(int i=0;i<lim;i+=(d<<1))
fo(j,0,d-1){
int tmp=a[i+j+d];
a[i+j+d]=(a[i+j]-tmp+mod)*tp%mod;
a[i+j]=(a[i+j]+tmp)*tp%mod;
}
}
高维前缀和的求法,枚举每一维,做一遍前缀和
在位运算上也就是枚举每一个二进制位,将这一位为\(0\)的加到这一位为\(1\)的上面
这样做高维前缀和虽然朴素,但是有正确性的保证,保证了不重不漏
code
#include<bits/stdc++.h>
using namespace std;
#define int long long
#define fo(i,x,y) for(int i=(x);i<=(y);i++)
#define fu(i,x,y) for(int i=(x);i>=(y);i--)
int read(){
int s=0,t=1;char ch=getchar();
while(!isdigit(ch)){if(ch=='-')t=-1;ch=getchar();}
while(isdigit(ch)){s=s*10+ch-'0';ch=getchar();}
return s*t;
}
const int mod=998244353;
const int N=1<<18;
int ksm(int x,int y){
int ret=1;
while(y){
if(y&1)ret=ret*x%mod;
x=x*x%mod;y>>=1;
}return ret;
}
int n,lim,a[N],b[N],c[N],aa[N],bb[N];
void init(int lim){fo(i,0,lim-1)a[i]=aa[i],b[i]=bb[i];}
void ort(int *a,int lim,int tp){
for(int d=1;d<lim;d<<=1)
for(int i=0;i<lim;i+=(d<<1))
fo(j,0,d-1){
a[i+j+d]=(a[i+j+d]+a[i+j]*tp)%mod;
}
}
void andt(int *a,int lim,int tp){
for(int d=1;d<lim;d<<=1)
for(int i=0;i<lim;i+=(d<<1))
fo(j,0,d-1){
a[i+j]=(a[i+j]+a[i+j+d]*tp)%mod;
}
}
void xort(int *a,int lim,int tp){
for(int d=1;d<lim;d<<=1)
for(int i=0;i<lim;i+=(d<<1))
fo(j,0,d-1){
int tmp=a[i+j+d];
a[i+j+d]=(a[i+j]-tmp+mod)*tp%mod;
a[i+j]=(a[i+j]+tmp)*tp%mod;
}
}
void mix(int lim){fo(i,0,lim-1)c[i]=a[i]*b[i]%mod;}
void pt(int *a,int lim){fo(i,0,lim-1)printf("%lld ",a[i]);printf("\n");}
signed main(){
n=read();lim=1<<n;
fo(i,0,lim-1)aa[i]=read();
fo(i,0,lim-1)bb[i]=read();
init(lim);ort(a,lim,1);ort(b,lim,1);mix(lim);ort(c,lim,mod-1);pt(c,lim);
init(lim);andt(a,lim,1);andt(b,lim,1);mix(lim);andt(c,lim,mod-1);pt(c,lim);
init(lim);xort(a,lim,1);xort(b,lim,1);mix(lim);xort(c,lim,ksm(2,mod-2));pt(c,lim);
}
分治FFT/NTT/FWT
这个东东就是把FFT和CDQ套在一起了
先处理左区间,处理当前区间,然后是右区间
code
#include<bits/stdc++.h>
using namespace std;
#define int long long
#define fo(i,x,y) for(int i=(x);i<=(y);i++)
#define fu(i,x,y) for(int i=(x);i>=(y);i++)
int read(){
int s=0,t=1;char ch=getchar();
while(!isdigit(ch)){if(ch=='-')t=-1;ch=getchar();}
while(isdigit(ch)){s=s*10+ch-'0';ch=getchar();}
return s*t;
}
const int mod=998244353;
const int N=1<<20;
int ksm(int x,int y){
int ret=1;
while(y){
if(y&1)ret=ret*x%mod;
x=x*x%mod;y>>=1;
}return ret;
}
int af[N],lim,len,w[N];
void ntt(int *a,int lim){
fo(i,0,lim-1)if(af[i]>i)swap(a[i],a[af[i]]);
for(int t=lim>>1,d=1;d<lim;d<<=1,t>>=1)
for(int i=0;i<lim;i+=(d<<1))
fo(j,0,d-1){
int tmp=w[t*j]*a[i+j+d]%mod;
a[i+j+d]=(a[i+j]-tmp+mod)%mod;
a[i+j]=(a[i+j]+tmp)%mod;
}
}
int n,f[N],g[N],a[N],b[N],c[N];
void mul(int n,int m,int v){
for(lim=1,len=0;lim<n+m;lim<<=1,len++);
fo(i,0,lim-1)af[i]=(af[i>>1]>>1)|((i&1)<<(len-1));
w[0]=1;w[1]=ksm(3,(mod-1)/lim);
fo(i,2,lim-1)w[i]=w[i-1]*w[1]%mod;
ntt(a,lim);ntt(b,lim);
w[0]=1;w[1]=ksm(w[1],mod-2);
fo(i,2,lim-1)w[i]=w[i-1]*w[1]%mod;
fo(i,0,lim-1)c[i]=a[i]*b[i]%mod;
ntt(c,lim);int iv=ksm(lim,mod-2);
fo(i,0,lim-1)c[i]=c[i]*iv%mod,a[i]=b[i]=0;
}
void sol(int l,int r){
// cout<<l<<" "<<r<<endl;
if(l==r)return ;
int mid=l+r>>1,ln=mid-l+1,lm=r-l+1;
sol(l,mid);
// cout<<"S "<<l<<" "<<r<<endl;
fo(i,l,mid)a[i-l]=f[i];//cout<<a[i]<<" ";cout<<endl;
fo(i,l,r)b[i-l]=g[i-l];//cout<<b[i]<<" ";cout<<endl;
mul(ln,lm,lm);
fo(i,ln,lm-1)f[i+l]=(f[i+l]+c[i])%mod;
sol(mid+1,r);
}
signed main(){
n=read();f[0]=1;
fo(i,1,n-1)g[i]=read();
sol(0,n-1);
fo(i,0,n-1)printf("%lld ",f[i]);
}
多项式求逆
就是求某一个多项式的逆元\(mod\ x^n\)意义下
主要吧,这个板子还挺难背的
递推上来的证明过程
所以结论就是\(B=2B'-AB'^2\)(\(B'\)表示\(mod\ x^{\lceil\frac{n}{2}\rceil}\)意义下的逆元)
code
#include<bits/stdc++.h>
using namespace std;
#define int long long
#define fo(i,x,y) for(int i=(x);i<=(y);i++)
#define fu(i,x,y) for(int i=(x);i>=(y);i--)
int read(){
int s=0,t=1;char ch=getchar();
while(!isdigit(ch)){if(ch=='-')t=-1;ch=getchar();}
while(isdigit(ch)){s=s*10+ch-'0';ch=getchar();}
return s*t;
}
const int mod=998244353;
const int N=1<<18;
int ksm(int x,int y){
int ret=1;
while(y){
if(y&1)ret=ret*x%mod;
x=x*x%mod;y>>=1;
}return ret;
}
int af[N],g[N],lim,len;
void ntt(int *a,int lim){
fo(i,0,lim-1)if(af[i]>i)swap(a[i],a[af[i]]);
for(int t=lim>>1,d=1;d<lim;d<<=1,t>>=1)
for(int i=0;i<lim;i+=(d<<1))
fo(j,0,d-1){
int tmp=g[t*j]*a[i+j+d]%mod;
a[i+j+d]=(a[i+j]-tmp+mod)%mod;
a[i+j]=(a[i+j]+tmp)%mod;
}
}
int a[N],b[N],c[N];
void calc(int deg){//deg表示当前mod的是x^deg
if(deg==1){return b[0]=ksm(a[0],mod-2),void();}
calc((deg+1)>>1);
for(lim=1,len=0;lim<deg*2;lim<<=1,len++);
fo(i,0,lim-1)af[i]=(af[i>>1]>>1)|((i&1)<<(len-1));
g[0]=1;g[1]=ksm(3,(mod-1)/lim);
fo(i,2,lim-1)g[i]=g[i-1]*g[1]%mod;
fo(i,0,deg-1)c[i]=a[i];//注意这里只赋值前deg-1位,因为这次乘法只能用到这些位
//这里a有deg-1位,而b只有deg+1>>1位,因为b是要平方的,并且这里是求a的deg-1位的逆元
ntt(c,lim);ntt(b,lim);
g[0]=1;g[1]=ksm(g[1],mod-2);
fo(i,2,lim-1)g[i]=g[i-1]*g[1]%mod;
fo(i,0,lim-1)b[i]=(2-c[i]*b[i]%mod+mod)*b[i]%mod;
ntt(b,lim);int iv=ksm(lim,mod-2);
fo(i,0,lim-1)b[i]=b[i]*iv%mod;
fo(i,deg,lim-1)b[i]=0;//把后边的清空,防止影响下一次运算
fo(i,0,lim-1)c[i]=0;//清空
}
int n;
signed main(){
n=read();
fo(i,0,n-1)a[i]=read();
calc(n);
fo(i,0,n-1)printf("%lld ",b[i]);
}
多项式开根
这个要用到上面的多项式求逆
这个就是求一个多项式的根号
于是乎我们仍然是递归求解
\(B=\frac{A+B^2}{2B}\) 证明
code
#include<bits/stdc++.h>
using namespace std;
#define int long long
#define fo(i,x,y) for(int i=(x);i<=(y);i++)
#define fu(i,x,y) for(int i=(x);i>=(y);i--)
int read(){
int s=0,t=1;char ch=getchar();
while(!isdigit(ch)){if(ch=='-')t=-1;ch=getchar();}
while(isdigit(ch)){s=s*10+ch-'0';ch=getchar();}
return s*t;
}
const int mod=998244353;
const int N=1<<20;
int ksm(int x,int y){
int ret=1;
while(y){
if(y&1)ret=ret*x%mod;
x=x*x%mod;y>>=1;
}return ret;
}
int af[N],lim,len,g[N],iv2;
void ntt(int *a,int lim){
fo(i,0,lim-1)if(af[i]>i)swap(a[i],a[af[i]]);
for(int t=lim>>1,d=1;d<lim;d<<=1,t>>=1)
for(int i=0;i<lim;i+=(d<<1))
fo(j,0,d-1){
int tmp=g[t*j]*a[i+j+d]%mod;
a[i+j+d]=(a[i+j]-tmp+mod)%mod;
a[i+j]=(a[i+j]+tmp)%mod;
}
}
int a[N],b[N],c[N];
void inv(int *x,int *y,int n){
if(n==1){return y[0]=ksm(x[0],mod-2),void();}
inv(x,y,n+1>>1);
for(lim=1,len=0;lim<2*n;lim<<=1,len++);
fo(i,0,lim-1)af[i]=(af[i>>1]>>1)|((i&1)<<(len-1));
g[0]=1;g[1]=ksm(3,(mod-1)/lim);
fo(i,2,lim-1)g[i]=g[i-1]*g[1]%mod;
fo(i,0,n-1)a[i]=x[i];
ntt(a,lim);ntt(y,lim);
g[0]=1;g[1]=ksm(g[1],mod-2);
fo(i,2,lim-1)g[i]=g[i-1]*g[1]%mod;
fo(i,0,lim-1)y[i]=(2ll-a[i]*y[i]%mod+mod)*y[i]%mod;
ntt(y,lim);int iv=ksm(lim,mod-2);
fo(i,0,lim-1)y[i]=y[i]*iv%mod,a[i]=0;
fo(i,n,lim-1)y[i]=0;
}
void sqt(int *x,int *y,int n){
if(n==1){return y[0]=1,void();}
sqt(x,y,n+1>>1);inv(y,c,n);
for(lim=1,len=0;lim<2*n;lim<<=1,len++);
fo(i,0,lim-1)af[i]=(af[i>>1]>>1)|((i&1)<<(len-1));
g[0]=1;g[1]=ksm(3,(mod-1)/lim);
fo(i,2,lim-1)g[i]=g[i-1]*g[1]%mod;
fo(i,0,n-1)a[i]=x[i];
fo(i,0,lim-1)c[i]=c[i]*iv2%mod;
ntt(a,lim);ntt(y,lim);ntt(c,lim);
g[0]=1;g[1]=ksm(g[1],mod-2);
fo(i,2,lim-1)g[i]=g[i-1]*g[1]%mod;
fo(i,0,lim-1)y[i]=(a[i]+y[i]*y[i]%mod)*c[i]%mod;
ntt(y,lim);int iv=ksm(lim,mod-2);
fo(i,0,lim-1)y[i]=y[i]*iv%mod,a[i]=b[i]=c[i]=0;
fo(i,n,lim-1)y[i]=0;
}
int f[N],h[N];
signed main(){
int n=read();iv2=ksm(2,mod-2);
fo(i,0,n-1)f[i]=read();
sqt(f,h,n);
fo(i,0,n-1)printf("%lld ",h[i]);
}
多项式ln
求导然后积分回去,完事
code
#include<bits/stdc++.h>
using namespace std;
#define int long long
#define fo(i,x,y) for(int i=(x);i<=(y);i++)
#define fu(i,x,y) for(int i=(x);i>=(y);i--)
int read(){
int s=0,t=1;char ch=getchar();
while(!isdigit(ch)){if(ch=='-')t=-1;ch=getchar();}
while(isdigit(ch)){s=s*10+ch-'0';ch=getchar();}
return s*t;
}
const int mod=998244353;
const int N=1<<20;
int ksm(int x,int y){
int ret=1;
while(y){
if(y&1)ret=ret*x%mod;
x=x*x%mod;y>>=1;
}return ret;
}
int af[N],lim,len,g[N];
void ntt(int *a,int lim){
fo(i,0,lim-1)if(af[i]>i)swap(a[i],a[af[i]]);
for(int t=lim>>1,d=1;d<lim;d<<=1,t>>=1)
for(int i=0;i<lim;i+=(d<<1))
fo(j,0,d-1){
int tmp=g[t*j]*a[i+j+d]%mod;
a[i+j+d]=(a[i+j]-tmp+mod)%mod;
a[i+j]=(a[i+j]+tmp)%mod;
}
}
int a[N],b[N],c[N];
void inv(int *x,int *y,int n){
if(n==1)return y[0]=ksm(x[0],mod-2),void();
inv(x,y,n+1>>1);
for(lim=1,len=0;lim<n*2;lim<<=1,len++);
fo(i,0,lim-1)af[i]=(af[i>>1]>>1)|((i&1)<<(len-1));
g[0]=1;g[1]=ksm(3,(mod-1)/lim);
fo(i,2,lim-1)g[i]=g[i-1]*g[1]%mod;
fo(i,0,n-1)a[i]=x[i];
ntt(a,lim);ntt(y,lim);
g[0]=1;g[1]=ksm(g[1],mod-2);
fo(i,2,lim-1)g[i]=g[i-1]*g[1]%mod;
fo(i,0,lim-1)y[i]=(2ll-a[i]*y[i]%mod+mod)*y[i]%mod;
ntt(y,lim);int iv=ksm(lim,mod-2);
fo(i,0,lim-1)y[i]=y[i]*iv%mod,a[i]=0;
fo(i,n,lim-1)y[i]=0;
}
int f[N],h[N],cf,ch;
void dec(int *x,int n){
fo(i,0,n-1)x[i]=x[i+1]*(i+1)%mod;
}
void sum(int *x,int n){
fu(i,n,1)x[i]=x[i-1]*ksm(i,mod-2)%mod;x[0]=0;
}
signed main(){
cf=read();
fo(i,0,cf-1)f[i]=read();
inv(f,h,cf);dec(f,cf);
for(lim=1,len=0;lim<=cf*2-1;lim<<=1,len++);
fo(i,0,lim-1)af[i]=(af[i>>1]>>1)|((i&1)<<(len-1));
g[0]=1;g[1]=ksm(3,(mod-1)/lim);
fo(i,2,lim-1)g[i]=g[i-1]*g[1]%mod;
ntt(h,lim);ntt(f,lim);
g[0]=1;g[1]=ksm(g[1],mod-2);
fo(i,2,lim-1)g[i]=g[i-1]*g[1]%mod;
fo(i,0,lim-1)h[i]=h[i]*f[i]%mod;
ntt(h,lim);int iv=ksm(lim,mod-2);
fo(i,0,lim-1)h[i]=h[i]*iv%mod;
sum(h,cf);fo(i,0,cf-1)printf("%lld ",h[i]);
}
多项式exp
牛顿迭代
code
#include<bits/stdc++.h>
using namespace std;
#define int long long
#define fo(i,x,y) for(int i=(x);i<=(y);i++)
#define fu(i,x,y) for(int i=(x);i>=(y);i--)
int read(){
int s=0,t=1;char ch=getchar();
while(!isdigit(ch)){if(ch=='-')t=-1;ch=getchar();}
while(isdigit(ch)){s=s*10+ch-'0';ch=getchar();}
return s*t;
}
const int mod=998244353;
const int N=1<<20;
int ksm(int x,int y){
int ret=1;
while(y){
if(y&1)ret=ret*x%mod;
x=x*x%mod;y>>=1;
}return ret;
}
int af[N],lim,len,w[N];
void ntt(int *a,int lim){
fo(i,0,lim-1)if(af[i]>i)swap(a[i],a[af[i]]);
for(int t=lim>>1,d=1;d<lim;d<<=1,t>>=1)
for(int i=0;i<lim;i+=(d<<1))
fo(j,0,d-1){
int tmp=w[t*j]*a[i+j+d]%mod;
a[i+j+d]=(a[i+j]-tmp+mod)%mod;
a[i+j]=(a[i+j]+tmp)%mod;
}
}
int a[N],b[N],c[N],d[N],e[N];
void inv(int *x,int *y,int n){
if(n==1)return y[0]=ksm(x[0],mod-2),void();
inv(x,y,n+1>>1);
for(lim=1,len=0;lim<n*2;lim<<=1,len++);
fo(i,0,lim-1)af[i]=(af[i>>1]>>1)|((i&1)<<(len-1));
w[0]=1;w[1]=ksm(3,(mod-1)/lim);
fo(i,2,lim-1)w[i]=w[i-1]*w[1]%mod;
fo(i,0,n-1)a[i]=x[i];
ntt(a,lim);ntt(y,lim);
w[0]=1;w[1]=ksm(w[1],mod-2);
fo(i,2,lim-1)w[i]=w[i-1]*w[1]%mod;
fo(i,0,lim-1)y[i]=(2ll-a[i]*y[i]%mod+mod)*y[i]%mod;
ntt(y,lim);int iv=ksm(lim,mod-2);
fo(i,0,lim-1)y[i]=y[i]*iv%mod,a[i]=0;
fo(i,n,lim-1)y[i]=0;
}
void dec(int *a,int n){
fo(i,0,n-2)a[i]=a[i+1]*(i+1)%mod;
a[n-1]=0;
}
void sum(int *a,int n){
fu(i,n-1,1)a[i]=a[i-1]*ksm(i,mod-2)%mod;
a[0]=0;
}
void ln(int *x,int *y,int n){
fo(i,0,n-1)c[i]=x[i];
dec(c,n);inv(x,b,n);
for(lim=1,len=0;lim<n*2;lim<<=1,len++);
fo(i,0,lim-1)af[i]=(af[i>>1]>>1)|((i&1)<<(len-1));
w[0]=1;w[1]=ksm(3,(mod-1)/lim);
fo(i,2,lim-1)w[i]=w[i-1]*w[1]%mod;
ntt(c,lim);ntt(b,lim);
w[0]=1;w[1]=ksm(w[1],mod-2);
fo(i,2,lim-1)w[i]=w[i-1]*w[1]%mod;
fo(i,0,lim-1)y[i]=c[i]*b[i]%mod;
ntt(y,lim);int iv=ksm(lim,mod-2);
fo(i,0,lim-1)y[i]=y[i]*iv%mod,c[i]=b[i]=0;
sum(y,n);fo(i,n,lim-1)y[i]=0;
}
void exp(int *x,int *y,int n){
if(n==1)return y[0]=1,void();
exp(x,y,n+1>>1);ln(y,d,n);
fo(i,0,n-1)e[i]=x[i];
for(lim=1,len=0;lim<n*2;lim<<=1,len++);
fo(i,0,lim-1)af[i]=(af[i>>1]>>1)|((i&1)<<(len-1));
w[0]=1;w[1]=ksm(3,(mod-1)/lim);
fo(i,2,lim-1)w[i]=w[i-1]*w[1]%mod;
ntt(e,lim);ntt(d,lim);ntt(y,lim);
w[0]=1;w[1]=ksm(w[1],mod-2);
fo(i,2,lim-1)w[i]=w[i-1]*w[1]%mod;
fo(i,0,lim-1)y[i]=(1ll-d[i]+e[i]+2*mod)*y[i]%mod;
ntt(y,lim);int iv=ksm(lim,mod-2);
fo(i,0,lim-1)y[i]=y[i]*iv%mod,e[i]=d[i]=0;
fo(i,n,lim-1)y[i]=0;
// cout<<n<<" ";fo(i,0,n-1)cout<<y[i]<<" ";cout<<endl;
}
int n,f[N],g[N];
signed main(){
n=read();
fo(i,0,n-1)f[i]=read();
exp(f,g,n);
fo(i,0,n-1)printf("%lld ",g[i]);
}
多项式除法
翻转多项式之后可以取模干掉余数,于是可以直接求逆然后带回去
code
#include<bits/stdc++.h>
using namespace std;
#define int long long
#define fo(i,x,y) for(int i=(x);i<=(y);i++)
#define fu(i,x,y) for(int i=(x);i>=(y);i--)
int read(){
int s=0,t=1;char ch=getchar();
while(!isdigit(ch)){if(ch=='-')t=-1;ch=getchar();}
while(isdigit(ch)){s=s*10+ch-'0';ch=getchar();}
return s*t;
}
const int mod=998244353;
const int N=1<<19;
int ksm(int x,int y){
int ret=1;
while(y){
if(y&1)ret=ret*x%mod;
x=x*x%mod;y>>=1;
}return ret;
}
int af[N],lim,len,w[N];
void ntt(int *a,int lim){
fo(i,0,lim-1)if(af[i]>i)swap(a[af[i]],a[i]);
for(int t=lim>>1,d=1;d<lim;d<<=1,t>>=1)
for(int i=0;i<lim;i+=(d<<1))
fo(j,0,d-1){
int tmp=w[t*j]*a[i+j+d]%mod;
a[i+j+d]=(a[i+j]-tmp+mod)%mod;
a[i+j]=(a[i+j]+tmp)%mod;
}
}
int a[N],b[N],c[N];
void inv(int *x,int *y,int n){
if(n==1)return y[0]=ksm(x[0],mod-2),void();
inv(x,y,n+1>>1);
for(lim=1,len=0;lim<n*2;lim<<=1,len++);
fo(i,0,lim-1)af[i]=(af[i>>1]>>1)|((i&1)<<(len-1));
w[0]=1;w[1]=ksm(3,(mod-1)/lim);
fo(i,2,lim-1)w[i]=w[i-1]*w[1]%mod;
fo(i,0,n-1)a[i]=x[i];
ntt(a,lim);ntt(y,lim);
w[0]=1;w[1]=ksm(w[1],mod-2);
fo(i,2,lim-1)w[i]=w[i-1]*w[1]%mod;
fo(i,0,lim-1)y[i]=(2ll-a[i]*y[i]%mod+mod)*y[i]%mod;
ntt(y,lim);int iv=ksm(lim,mod-2);
fo(i,0,lim-1)y[i]=y[i]*iv%mod,a[i]=0;
fo(i,n,lim-1)y[i]=0;
}
void mul(int *x,int *y,int n,int m,int *z,int v){
for(lim=1,len=0;lim<n+m;lim<<=1,len++);
fo(i,0,lim-1)af[i]=(af[i>>1]>>1)|((i&1)<<(len-1));
w[0]=1;w[1]=ksm(3,(mod-1)/lim);
fo(i,2,lim-1)w[i]=w[i-1]*w[1]%mod;
fo(i,0,n-1)a[i]=x[i];
fo(i,0,m-1)b[i]=y[i];
ntt(a,lim);ntt(b,lim);
w[0]=1;w[1]=ksm(w[1],mod-2);
fo(i,2,lim-1)w[i]=w[i-1]*w[1]%mod;
fo(i,0,lim-1)z[i]=a[i]*b[i]%mod;
ntt(z,lim);int iv=ksm(lim,mod-2);
fo(i,0,lim-1)z[i]=z[i]*iv%mod,a[i]=b[i]=0;
fo(i,v,lim-1)z[i]=0;
}
int n,m,f[N],g[N],ng[N],q[N],r[N];
signed main(){
n=read();m=read();
fo(i,0,n)f[i]=read();
fo(i,0,m)g[i]=read();
fo(i,0,n>>1)swap(f[i],f[n-i]);
fo(i,0,m>>1)swap(g[i],g[m-i]);
inv(g,ng,n-m+1);
mul(f,ng,n+1,n-m+1,q,n-m+1);
fo(i,0,n-m>>1)swap(q[i],q[n-m-i]);
fo(i,0,n>>1)swap(f[i],f[n-i]);
fo(i,0,m>>1)swap(g[i],g[m-i]);
mul(q,g,n-m+1,m+1,ng,m);
fo(i,0,m-1)r[i]=(f[i]-ng[i]+mod)%mod;
fo(i,0,n-m)printf("%lld ",q[i]);printf("\n");
fo(i,0,m-1)printf("%lld ",r[i]);
}
技巧
1、生成函数时,有关于对称的问题,如果对称的话,那么两个点的位置相加是对称中心的两倍,意思就是多项式卷积之后对称轴是同一条的都在一起
2、对于生成函数的利用,我们只需要系数,而并不是真正的把\(x\)带进去找值,我只需要系数就行了
比如:\(\sum_{i=0}^{n}\sum_{j=0}^{i}j!*2^{i-j}\),这个东西阶乘和幂可以分别看作是两个多项式的系数
这样我们可以通过卷积快速求出后面的和,也就是卷积后多项式的第\(i\)项系数,符合多项式乘法的规则,然后带入求和