多项式求逆/分治FFT 学习笔记

一、多项式求逆

  • 给定一个多项式 \(F(x)\),请求出一个多项式 \(G(x)\), 满足 \(F(x) * G(x) \equiv 1 ( \mathrm{mod\:} x^n )\)。系数对 \(998244353\)取模。
  • 考虑递归求解,当\(F\)的最高次为\(0\)时,\(G_0=F_0^{-1}\)
  • 假设我们知道了\(F(x)\)在模\(x^{\left \lceil \frac{n}{2}\right \rceil}\)意义下的逆元\(G'\)
  • 那么\(F∗G′≡1(\mathrm{mod\:} x^{\left \lceil \frac{n}{2}\right \rceil})\)\(F∗G≡1(\mathrm{mod\:} x^{\left \lceil \frac{n}{2}\right \rceil})\)
  • 因此 \(G-G'\equiv 0(\mathrm{mod\:} x^{\left \lceil \frac{n}{2}\right \rceil})\)
  • 然后两边平方:\((G-G')^2\equiv 0(\mathrm{mod\:} x^n)\)
  • 所以\(G^2-2GG'+G'^2\equiv0(\mathrm{mod\:} x^n)\)
  • 通乘\(F\)后,由于\(F*G \equiv0(\mathrm{mod\:} x^n)\),所以\(G-2G'+FG'^2\equiv0(\mathrm{mod\:} x^n)\)
  • 最后得到\(G\equiv2G'-FG'^2(\mathrm{mod\:} x^n)\)
  • 总复杂度\(O(n \log n)\)

二、分治FFT

  • 在求卷积的时候,如果后面的数字基于前面的数字,朴素的\(FFT\)就会退化至\(O(n^2 \log n)\)
  • 考虑\(cdq\)分治
  • 我们先求出\(l \rightarrow mid\)的答案,然后考虑它对$mid+1 \rightarrow r $的贡献。
  • 显然,对于\(mid+1 \rightarrow r\)\(x\)贡献是:

\[ \sum_{i=l}^{mid}f[i]g[x-i] \]

  • 贡献可以用\(FFT\)计算
  • 总复杂度\(O(n \log n)\)

传送门luoguP4721 【模板】分治FFT

Description

给定长度为 \(n-1\)的数组$ g[1],g[2],..,g[n-1]g[1],g[2],..,g[n−1]$,求 \(f[0],f[1],..,f[n-1]f[0],f[1],..,f[n−1]\),其中

\[f[i]=\sum_{j=1}^if[i-j]g[j]\]

边界为 \(f[0]=1\) 。答案模 \(998244353\)


Code-多项式求逆版 

由题可知:
\[ f*g=f-1 \]
所以:
\[ f\equiv(1-g)^{-1}(\mathrm{mod\:} x^n) \]
直接多项式求逆就可以啦

复杂度\(O(n\log n)\)

#include<bits/stdc++.h>
#define ll long long
#define max(a,b) ((a)>(b)?(a):(b))
#define min(a,b) ((a)<(b)?(a):(b))
#define swap(a,b) (a^=b^=a^=b)
inline ll read()
{
    ll x=0,f=1;char ch=getchar();
    while(ch<'0'||ch>'9'){if(ch=='-')f=-1;ch=getchar();}
    while(ch>='0'&&ch<='9'){x=(x<<3)+(x<<1)+ch-'0';ch=getchar();}
    return x*f;
}
#define mod 998244353
#define g 3
#define invg 332748118
#define MN 2097152
ll a[MN],b[MN],c[MN],N,di,invN;
int pos[MN];
bool now;
inline ll fpow(ll x,int m){ll res=1;for(;m;m>>=1,x=x*x%mod) if(m&1)res=res*x%mod;return res;}
inline void NTT(ll *a,int type)
{
    register int i,j,p,k;
    for(i=0;i<N;++i)if(i<pos[i]) swap(a[i],a[pos[i]]);
    for(i=1;i<N;i<<=1)
    {
        ll wn=fpow(type>0?g:invg,(mod-1)/(i<<1));
        for(p=i<<1,j=0;j<N;j+=p) 
        {
            ll w=1;
            for(k=0;k<i;++k,w=w*wn%mod)
            {
                ll X=a[j+k],Y=w*a[j+i+k]%mod;
                a[j+k]=(X+Y)%mod;a[j+i+k]=(X-Y+mod)%mod;
            }
        }
    }
}
void solve(int n,ll *a,ll *b)
{
    if(n==1){b[0]=fpow(a[0],mod-2);return;}
    solve((n+1)>>1,a,b);
    for(N=1,di=0;N<(n<<1);N<<=1,++di);
    register int i;invN=fpow(N,mod-2);
    for(i=0;i<N;++i) pos[i]=(pos[i>>1]>>1)|((i&1)<<(di-1));
    for(i=0;i<N;++i) c[i]=a[i]*(i<n);
    NTT(c,1),NTT(b,1);
    for(i=0;i<N;++i) b[i]=1ll*(2-1ll*c[i]*b[i]%mod+mod)%mod*b[i]%mod;
    NTT(b,-1);for(i=0;i<N;++i) b[i]=b[i]*invN%mod;
    for(i=n;i<N;++i) b[i]=0;
}
int main()
{
    register int n=read(),i;
    for(a[0]=1,i=1;i<n;++i) a[i]=(mod-read())%mod;
    solve(n,a,b);
    for(i=0;i<n;++i) printf("%lld ",b[i]);
    return 0;
}


Code-分治FFT版 

复杂度\(O(n\log^2 n)\)

#include<bits/stdc++.h>
#define ll long long
#define max(a,b) ((a)>(b)?(a):(b))
#define min(a,b) ((a)<(b)?(a):(b))
inline int read()
{
    int x=0,f=1;char ch=getchar();
    while(ch<'0'||ch>'9'){if(ch=='-')f=-1;ch=getchar();}
    while(ch>='0'&&ch<='9'){x=(x<<3)+(x<<1)+ch-'0';ch=getchar();}
    return x*f;
}
#define mod 998244353
#define g 3
#define invg 332748118
#define MN 262144
ll G[MN],F[MN],N,di,pos[MN],A[MN],B[MN],invN;
inline ll fpow(ll x,int m){ll res=1;for(;m;m>>=1,x=x*x%mod) if(m&1)res=res*x%mod;return res;}
inline void NTT(ll *a,int type)
{
    register int i,j,p,k;
    for(i=0;i<N;++i)if(i<pos[i]) std::swap(a[i],a[pos[i]]);
    for(i=1;i<N;i<<=1)
    {
        ll wn=fpow(type>0?g:invg,(mod-1)/(i<<1));
        for(p=i<<1,j=0;j<N;j+=p) 
        {
            ll w=1;
            for(k=0;k<i;++k,w=w*wn%mod)
            {
                ll X=a[j+k],Y=w*a[j+i+k]%mod;
                a[j+k]=(X+Y)%mod;a[j+i+k]=(X-Y+mod)%mod;
            }
        }
    }
    if(type==-1) for(i=0;i<N;++i) a[i]=a[i]*invN%mod;
}
inline void cdqNTT(ll *a,ll *b,int l,int r)
{
    if(l==r) return;
    register int mid=(l+r)>>1,i,len=r-l+1;
    cdqNTT(a,b,l,mid);
    
    for(N=1,di=0;N<len<<1;N<<=1,++di);
    for(i=0;i<N;++i) pos[i]=(pos[i>>1]>>1)|((i&1)<<(di-1)); 
    invN=fpow(N,mod-2);
    
    for(i=0;i<N;++i) A[i]=B[i]=0;
    for(i=l;i<=mid;++i) A[i-l]=a[i];
    for(i=0;i<=r-l;++i) B[i]=b[i];
    NTT(A,1),NTT(B,1);
    for(i=0;i<N;++i) A[i]=A[i]*B[i]%mod;
    NTT(A,-1);
    
    for(int i=mid+1;i<=r;++i) (a[i]+=A[i-l])%=mod;
    cdqNTT(a,b,mid+1,r);
}
int main()
{
    register int n,i;
    n=read();
    for(i=1;i<n;++i) G[i]=read();
    F[0]=1;
    cdqNTT(F,G,0,n-1);
    for(i=0;i<n;++i) printf("%lld ",F[i]);
    return 0;
}



Blog来自PaperCloud,未经允许,请勿转载,TKS!

转载于:https://www.cnblogs.com/PaperCloud/p/10254328.html

  • 0
    点赞
  • 0
    收藏
    觉得还不错? 一键收藏
  • 0
    评论
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值