多项式全家桶

分类: IT文章 • 2022-05-05 13:45:02

多项式乘法：

然而蒟蒻的我并不会证明

$FFT$：

 1 struct cp{
 2 dd x,y;
 3 friend cp operator + (const cp &s1,const cp &s2){ return (cp){s1.x+s2.x,s1.y+s2.y}; }
 4 friend cp operator - (const cp &s1,const cp &s2){ return (cp){s1.x-s2.x,s1.y-s2.y}; }
 5 friend cp operator * (const cp &s1,const cp &s2){ return (cp){s1.x*s2.x-s1.y*s2.y,s1.y*s2.x+s1.x*s2.y}; }
 6 }A[N1],B[N1],C[N1];
 7  
 8 int r[N1];
 9 void FFT(cp *s,int len,int type)
10 {
11     int i,j,k;
12     for(i=0;i<len;i++) if(i<r[i]) swap(s[i],s[r[i]]);
13     for(k=2;k<=len;k<<=1)
14     {
15         cp wn=(cp){cos(2.0*pi*type/k),sin(2.0*pi*type/k)},w,t;
16         for(i=0;i<len;i+=k)
17         {
18             w=(cp){1,0};
19             for(j=0;j<(k>>1);j++,w=w*wn)
20             {
21                 t=w*s[i+j+(k>>1)];
22                 s[i+j+(k>>1)]=s[i+j]-t;
23                 s[i+j]=s[i+j]+t;
24             }
25         }
26     }
27 }
28 void FFT_Main(int len)
29 {
30     FFT(A,len,1); FFT(B,len,1);
31     for(int i=0;i<len;i++) C[i]=A[i]*B[i];
32     FFT(C,len,-1);
33     for(int i=0;i<len;i++) C[i].x/=len;
34 }
35

View Code

$NTT$：

 1 ll qpow(ll x,ll y,const ll &mod)
 2 {
 3     ll ans=1;
 4     for(;y;x=x*x%mod,y>>=1) 
 5         if(y&1) ans=ans*x%mod;
 6     return ans;
 7 }
 8 const ll p=998244353;
 9 int r[N1];
10 ll A[N1],B[N1],C[N1],mulwn[N1],invwn[N1],invl,inv2;
11 void Pre(int len,int L)
12 {
13     int i; ll inv,invy;
14     for(i=0;i<len;i++) r[i]=(r[i>>1]>>1)|((i&1)<<(L-1));
15     for(i=1;i<=len;i<<=1) 
16     {
17         mulwn[i]=qpow(3,(p-1)/i,p);
18         exgcd(mulwn[i],p,inv,invy);
19         invwn[i]=(inv%p+p)%p;
20     }
21     exgcd(len,p,invl,invy); invl=(invl%p+p)%p;
22 }
23 void NTT(ll *s,int len,int type)
24 {
25     int i,j,k; ll w,t,wn;
26     for(i=0;i<len;i++) if(i<r[i]) swap(s[i],s[r[i]]);
27     for(k=2;k<=len;k<<=1)
28     {
29         wn=type>0?mulwn[k]:invwn[k];
30         for(i=0;i<len;i+=k)
31         {
32             for(j=0,w=1;j<(k>>1);j++,w=w*wn%p)
33             {
34                 t=w*s[i+j+(k>>1)]%p;
35                 s[i+j+(k>>1)]=(s[i+j]-t+p)%p;
36                 s[i+j]=(s[i+j]+t)%p;
37             }
38         }
39     }
40 }
41 void NTT_Main(int len)
42 {
43     NTT(A,len,1); NTT(B,len,1);
44     for(int i=0;i<len;i++) C[i]=A[i]*B[i]%p;
45     NTT(C,len,-1);
46     for(int i=0;i<len;i++) C[i]=C[i]*invl%p;
47 }

View Code

多项式$+CDQ$分治：

有一些卷积形式，后面的项需要依赖前面项的答案

比如$f(i)=sum_{j=0}^{i}g(j)f(i-j)$

所以我们用$CDQ$分治处理前面对后面项的影响

即用$[l,mid]$的答案更新$(mid,r]$

分治$NTT$：

 1 const ll p=998244353;
 2 ll qpow(ll x,ll y,const ll &mod)
 3 {
 4     ll ans=1;
 5     while(y){
 6         if(y&1) ans=ans*x%mod;
 7         x=x*x%mod; y>>=1;
 8     }return ans;
 9 }
10 int r[19][N1];
11 ll A[N1],B[N1],C[N1],g[N1],f[N1],ret[N1],invl[N1],mulwn[N1],invwn[N1];
12 void Pre(int len,int L)
13 {
14     int i,j;ll inv,invy;
15     for(j=1;j<=L;j++) for(i=0;i<(1<<j);i++) 
16         r[j][i]=(r[j][i>>1]>>1)|((i&1)<<(j-1)); 
17     for(i=2;i<=len;i<<=1) 
18     {
19         mulwn[i]=qpow(3,(p-1)/i,p); 
20         exgcd(mulwn[i],p,inv,invy); invwn[i]=(inv%p+p)%p;
21         exgcd(i,p,inv,invy); invl[i]=(inv%p+p)%p;
22     }
23 }
24 void NTT(ll *s,int len,int type,int L)
25 {
26     int i,j,k,inv; ll w,wn,t;
27     for(i=0;i<len;i++) 
28         if(i<r[L][i]) swap(s[i],s[r[L][i]]);
29     for(k=2;k<=len;k<<=1)
30     {
31         wn=type>0?mulwn[k]:invwn[k];
32         for(i=0;i<len;i+=k)
33         {
34             for(j=0,w=1;j<(k>>1);j++,w=w*wn%p)
35             {
36                 t=w*s[i+j+(k>>1)]%p;
37                 s[i+j+(k>>1)]=(s[i+j]-t+p)%p;
38                 s[i+j]=(s[i+j]+t)%p;
39             }
40         }
41     }
42     if(type==-1)
43     for(i=0;i<len;i++)
44         s[i]=s[i]*invl[len]%p;
45 }
46 void NTT_Main(int len,int L)
47 {
48     NTT(A,len,1,L); NTT(B,len,1,L);
49     for(int i=0;i<len;i++) C[i]=A[i]*B[i]%p;
50     NTT(C,len,-1,L);
51 }
52 void CDQ(int l,int r,int L)
53 {
54     if(l==r){ ret[l]=f[l]; return; }
55     int mid=(l+r)>>1,i;
56     CDQ(l,mid,L-1);
57     for(i=l;i<=r;i++) A[i-l]=ret[i],B[i-l]=g[i-l];
58     NTT_Main(r-l+1,L);
59     for(i=mid+1;i<=r;i++) (f[i]+=C[i-l])%=p;
60     CDQ(mid+1,r,L-1);
61 }
62 
63 int n,m,len,L;
64 int main()
65 {
66     scanf("%d",&n); int i; f[0]=1;
67     for(i=1;i<n;i++) g[i]=gint(); 
68     for(len=1,L=0;len<n;len<<=1,L++);
69     Pre(len,L);
70     CDQ(0,len-1,L);
71     for(i=0;i<n;i++) printf("%lld ",f[i]);
72     puts("");
73     return 0;
74 }

View Code

多项式求逆：

$mod;x^{n}$表示取出这个多项式的0~n-1次项

假设我们已经求出了$H(x)$，满足$H(x)F(x)equiv 1;(mod;x^{ left lceil frac{n}{2} ight ceil })$

我们想求$G(x)F(x)equiv 1;(mod;x^{n}) $

$(G(x)-H(x))equiv 0;(mod;x^{ left lceil frac{n}{2} ight ceil })$

$(G(x)-H(x))^{2} equiv 0;(mod;x^{n})$

展开

$G^{2}(x)-2G(x)H(x)+H^{2}(x) equiv 0;(mod;x^{n})$

$G(x)=2H(x)-frac{H^{2}(x)}{G(x)}$

$=2H(x)-H^{2}(x)F(x)$

$=H(x)*(2-H(x)F(x));(mod;x^{n})$

然后上$NTT$就好了

 1 void NTT(ll *s,int len,int type)
 2 {
 3     int i,j,k; ll wn,w,t;
 4     for(i=0;i<len;i++) if(i<r[i]) swap(s[i],s[r[i]]);
 5     for(k=2;k<=len;k<<=1)
 6     {
 7         wn=type>0?mulwn[k]:invwn[k];
 8         for(i=0;i<len;i+=k)
 9         {
10             for(j=0,w=1;j<(k>>1);j++,w=w*wn%p)
11             {
12                 t=w*s[i+j+(k>>1)]%p;
13                 s[i+j+(k>>1)]=(s[i+j]+p-t)%p;
14                 s[i+j]=(s[i+j]+t)%p;
15             }
16         }
17     }
18     if(type==-1)
19     {
20         ll invl=qpow(len,p-2);
21         for(i=0;i<len;i++) s[i]=s[i]*invl%p;
22     }
23 }
24 int len,L;
25 void NTT_INV(int n)
26 {
27     if(n==1){ b[0]=qpow(f[0],p-2); return; }
28     NTT_INV((n+1)>>1); int i;
29     for(len=1,L=0;len<n+2*((n+1)>>1);len<<=1,L++);
30     for(i=0;i<len;i++) r[i]=(r[i>>1]>>1)|((i&1)<<(L-1));
31     memset(a,0,len<<3); memcpy(a,f,n<<3);
32     NTT(a,len,1); NTT(b,len,1);
33     for(i=0;i<len;i++)
34         b[i]=(2ll-a[i]*b[i]%p+p)%p*b[i]%p;
35     NTT(b,len,-1);
36     for(i=n;i<len;i++) b[i]=0;
37 }
38 
39 int n;
40 
41 int main()
42 {
43     int ret=0,i; 
44     scanf("%d",&n);
45     for(i=0;i<n;i++) scanf("%lld",&f[i]);
46     for(len=1,L=0;len<n+n-1;len<<=1,L++);
47     for(i=1;i<=len;i<<=1) mulwn[i]=qpow(3,(p-1)/i), invwn[i]=qpow(mulwn[i],p-2);
48     NTT_INV(n);
49     for(i=0;i<n;i++) printf("%lld ",b[i]);
50     puts("");
51     return 0;
52 }

View Code

多项式求ln：

根据求导公式，若$f(x)=ln;x$，那么$f'(x)=frac{1}{x}$

现在已知多项式$F(x)$，那么$ln;F(x)$的一阶导$=frac{1}{F'(x)}$

我们利用导数知识求出$F'(x)$，$x^{a}$的一阶导是$ax^{a-1}$，$F'(x)$就是对$x$的不同次项求导，再相加

再用多项式求逆，求出多项式$G(x)=frac{1}{F'(x)}$，再通过不定积分推回去

此处的不定积分其实就是把导数公式$x^{a}=ax^{a-1}$逆推回去

 1 ll qpow(ll x,ll y)
 2 {
 3     ll ans=1;
 4     for(;y;x=x*x%p,y>>=1) 
 5         if(y&1) ans=ans*x%p;
 6     return ans;
 7 }
 8 
 9 int r[N1];
10 ll a[N1],b[N1],c[N1],mulwn[N1],invwn[N1],f[N1];
11 
12 void NTT(ll *s,int len,int type)
13 {
14     int i,j,k; ll wn,w,t;
15     for(i=0;i<len;i++) if(i<r[i]) swap(s[i],s[r[i]]);
16     for(k=2;k<=len;k<<=1)
17     {
18         wn=type>0?mulwn[k]:invwn[k];
19         for(i=0;i<len;i+=k)
20         {
21             for(j=0,w=1;j<(k>>1);j++,w=w*wn%p)
22             {
23                 t=w*s[i+j+(k>>1)]%p;
24                 s[i+j+(k>>1)]=(s[i+j]+p-t)%p;
25                 s[i+j]=(s[i+j]+t)%p;
26             }
27         }
28     }
29     if(type==-1)
30     {
31         ll invl=qpow(len,p-2);
32         for(i=0;i<len;i++) s[i]=s[i]*invl%p;
33     }
34 }
35 int len,L;
36 void NTT_INV(int n)
37 {
38     if(n==1){ b[0]=qpow(f[0],p-2); return; }
39     NTT_INV((n+1)>>1); int i;
40     for(len=1,L=0;len<n+2*((n+1)>>1);len<<=1,L++);
41     for(i=0;i<len;i++) r[i]=(r[i>>1]>>1)|((i&1)<<(L-1));
42     memset(a,0,len<<3); memcpy(a,f,n<<3);
43     NTT(a,len,1); NTT(b,len,1);
44     for(i=0;i<len;i++)
45         b[i]=(2ll-a[i]*b[i]%p+p)%p*b[i]%p;
46     NTT(b,len,-1);
47     for(i=n;i<len;i++) b[i]=0;
48 }
49 
50 int n;
51 
52 int main()
53 {
54     int ret=0,i; 
55     scanf("%d",&n);
56     for(i=0;i<n;i++) scanf("%lld",&f[i]);
57     for(len=1,L=0;len<n+n-1;len<<=1,L++);
58     for(i=1;i<=len;i<<=1) mulwn[i]=qpow(3,(p-1)/i), invwn[i]=qpow(mulwn[i],p-2);
59     NTT_INV(n);
60     for(i=0;i<n;i++) printf("%lld ",b[i]);
61     puts("");
62     return 0;
63 }

View Code

多项式全家桶

相关推荐