人傻常数大.jpg
因为求逆的时候没清零结果调了几个小时……
前置芝士
多项式除法,多项式求逆
什么?你不会?左转你谷模板区,包教包会
题解
首先我们要知道一个结论\[f(x_0)\equiv f(x)\pmod{(x-x_0)}\]
其中\(x_0\)为一个常量,\(f(x_0)\)也为一个常量
证明如下,设\(f(x)=g(x)(x-x_0)+A\),也就是说\(A\)是\(f(x)\)对\((x-x_0)\)这个多项式取模之后的结果
因为\((x-x_0)\)的最高次项为\(1\),所以\(A\)的最高次项为\(0\),也就是说\(A\)是一个常数,即\(f(x)\equiv A\pmod{(x-x_0)}\)
我们把\(x_0\)代入上式,得\(f(x_0)=g(x_0)(x_0-x_0)+A\),同理可得\(f(x_0)\equiv A\pmod{(x-x_0)}\)
于是我们知道上式成立
这有毛用啊\(O(n\log n)\)多项式取模还没我暴力快
乍一看的确没啥卵用,但是考虑取模的过程是否能优化呢?
答案是可以的,我们考虑分治。设当前分治区间为\([l,r]\),令\(P_0=\prod_{i=l}^{mid}(x-x_i)\),\(P_1=\prod_{i=mid+1}^r (x-x_i)\),当前已经算出了\(A\equiv f(x)\pmod{\prod_{i=l}^r(x-x_i)}\),那么只要分别用\(A\)对\(P_0\)和\(P_1\)取模,然后继续递归下去就行了。取模之后\(A(x)\)的最高次项的次数变为原来的一半,问题规模也就变为原来的一半。继续递归下去就行了
时间复杂度为\(O(n\log^2n)\)
upd:改了改代码,常数应该会小一点,比方说分治到某个时候暴力秦九韶展开
//minamoto#include#define R register#define fp(i,a,b) for(R int i=(a),I=(b)+1;i I;--i)#define go(u) for(int i=head[u],v=e[i].v;i;i=e[i].nx,v=e[i].v)using namespace std;char buf[1<<21],*p1=buf,*p2=buf;inline char getc(){return p1==p2&&(p2=(p1=buf)+fread(buf,1,1<<21,stdin),p1==p2)?EOF:*p1++;}int read(){ R int res,f=1;R char ch; while((ch=getc())>'9'||ch<'0')(ch=='-')&&(f=-1); for(res=ch-'0';(ch=getc())>='0'&&ch<='9';res=res*10+ch-'0'); return res*f;}char sr[1<<21],z[20];int C=-1,Z=0;inline void Ot(){fwrite(sr,1,C+1,stdout),C=-1;}void print(R int x){ if(C>1<<20)Ot();if(x<0)sr[++C]='-',x=-x; while(z[++Z]=x%10+48,x/=10); while(sr[++C]=z[Z],--Z);sr[++C]='\n';}const int N=(1<<17)+5,P=998244353;inline int add(R int x,R int y){return x+y>=P?x+y-P:x+y;}inline int dec(R int x,R int y){return x-y<0?x-y+P:x-y;}inline int mul(R int x,R int y){return 1ll*x*y-1ll*x*y/P*P;}int ksm(R int x,R int y){ R int res=1; for(;y;y>>=1,x=mul(x,x))(y&1)?res=mul(res,x):0; return res;}int r[19][N],w[2][N],lg[N],inv[19];void Pre(){ fp(d,1,17){ fp(i,1,(1< >1]>>1)|((i&1)<<(d-1)); lg[1< < >1,i=1,x,y;i<131072;i<<=1,t>>=1){ x=ksm(3,t),y=ksm(332748118,t),w[0][i]=w[1][i]=1; fp(k,1,i-1) w[1][k+i]=mul(w[1][k+i-1],x), w[0][k+i]=mul(w[0][k+i-1],y); }}int lim,d,n,m;inline void init(R int len){lim=1,d=0;while(lim <<=1,++d;}void NTT(int *A,int ty){ fp(i,0,lim-1)if(i >1),lim=(len<<1),d=lg[lim]; static int A[N],B[N]; fp(i,0,len-1)A[i]=a[i],B[i]=b[i];fp(i,len,lim-1)A[i]=B[i]=0; NTT(A,1),NTT(B,1); fp(i,0,lim-1)A[i]=mul(A[i],mul(B[i],B[i])); NTT(A,0); fp(i,0,len-1)b[i]=dec(add(b[i],b[i]),A[i]); fp(i,len,lim-1)b[i]=0;}struct node{ node *lc,*rc;vector vec;int deg; void Mod(const int *a,int *r,int n){ static int A[N],B[N],D[N]; int len=1;while(len<=n-deg)len<<=1; fp(i,0,n)A[i]=a[n-i];fp(i,0,deg)B[i]=vec[deg-i]; fp(i,n-deg+1,len-1)B[i]=0; Inv(B,D,len); lim=(len<<1),d=lg[lim]; fp(i,n-deg+1,lim-1)A[i]=D[i]=0; NTT(A,1),NTT(D,1); fp(i,0,lim-1)A[i]=mul(A[i],D[i]); NTT(A,0); reverse(A,A+n-deg+1); init(n+1); fp(i,n-deg+1,lim-1)A[i]=0; fp(i,0,deg)B[i]=vec[i];fp(i,deg+1,lim-1)B[i]=0; NTT(A,1),NTT(B,1); fp(i,0,lim-1)A[i]=mul(A[i],B[i]); NTT(A,0); fp(i,0,deg-1)r[i]=dec(a[i],A[i]); } void Mul(){ static int A[N],B[N];deg=lc->deg+rc->deg,vec.resize(deg+1),init(deg+1); fp(i,0,lc->deg)A[i]=lc->vec[i];fp(i,lc->deg+1,lim-1)A[i]=0; fp(i,0,rc->deg)B[i]=rc->vec[i];fp(i,rc->deg+1,lim-1)B[i]=0; NTT(A,1),NTT(B,1); fp(i,0,lim-1)A[i]=mul(A[i],B[i]); NTT(A,0); fp(i,0,deg)vec[i]=A[i]; }}pool[N],*rt;int A[N],a[N],tot;inline node* newnode(){return &pool[tot++];}void solve(node* &p,int l,int r){ p=newnode(); if(l==r)return p->deg=1,p->vec.resize(2),p->vec[0]=P-a[l],p->vec[1]=1,void(); int mid=(l+r)>>1; solve(p->lc,l,mid),solve(p->rc,mid+1,r); p->Mul();}int b[25];void calc(node* p,int l,int r,const int *A){ if(r-l<=512){ fp(i,l,r){ int x=a[i],c1,c2,c3,c4,now=A[r-l]; b[0]=1;fp(j,1,16)b[j]=mul(b[j-1],x); for(R int j=r-l-1;j-15>=0;j-=16){ c1=(1ll*now*b[16]+1ll*A[j]*b[15]+1ll*A[j-1]*b[14]+1ll*A[j-2]*b[13])%P, c2=(1ll*A[j-3]*b[12]+1ll*A[j-4]*b[11]+1ll*A[j-5]*b[10]+1ll*A[j-6]*b[9])%P, c3=(1ll*A[j-7]*b[8]+1ll*A[j-8]*b[7]+1ll*A[j-9]*b[6]+1ll*A[j-10]*b[5])%P, c4=(1ll*A[j-11]*b[4]+1ll*A[j-12]*b[3]+1ll*A[j-13]*b[2]+1ll*A[j-14]*b[1])%P, now=(0ll+c1+c2+c3+c4+A[j-15])%P; } fd(j,(r-l)%16-1,0)now=(1ll*now*x+A[j])%P; print(now); } return; } int mid=(l+r)>>1,b[p->deg+1]; p->lc->Mod(A,b,p->deg-1),calc(p->lc,l,mid,b); p->rc->Mod(A,b,p->deg-1),calc(p->rc,mid+1,r,b);}int main(){// freopen("testdata.in","r",stdin); n=read(),m=read();if(!m)return 0; Pre(); fp(i,0,n)A[i]=read(); fp(i,1,m)a[i]=read(); solve(rt,1,m); if(n>=m)rt->Mod(A,A,n); calc(rt,1,m,A); return Ot(),0;}