用户:
爵士OIer查看:6 回复:24 评论:6 创建时间:2021-04-10T14:27:13
本期内容比较简单而实用,欢迎学习!
一个月前我向大家分享了快速傅里叶变换和快速数论变换:
https://shequ.codemao.cn/community/3喵745
现在让我们继续多项式科技的进程!注:内容建立在上一期的基础上。
因此我非常良心地给大家每期教程都回顾一遍FFT和NTT!
点值表示法
众所周知,坐标系内 n 个不同的点确定 n−1 次的多项式。

单位根



原根

我们可以用原根替代单位根来实现NTT。
至于代码大家就去上面那个链接看啦qwq
多项式求逆
如果我们在数据运算的时候遇到包含了下面这个问题,我们怎么解决?

这个问题就是多项式求逆问题。求解它的主要思想是倍增,递归实现。
倍增过程

代码实现Code
基本上是很好写的。
时间复杂度 O (n log n)。

#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
const ll mod=998244353;
#define G 3
#define MAXN 400040
#define invG 332748118
#define ck(x) ((x)>=mod?(x)-mod:(x))
ll n,tr[MAXN],f[MAXN],g[MAXN],h[MAXN];
inline ll read()
{
int q=0;char ch=' ';
while(ch<'0'||ch>'9')ch=getchar();
while(ch>='0'&&ch<='9')q=(q<<3)+(q<<1)+ch-'0',ch=getchar();
return q;
}
inline ll Quickly_Power(ll a,ll b)
{
ll ans=1;
while(b)
{
if(b&1)ans=ans*a%mod;
a=a*a%mod,b>>=1;
}
return ans;
}
void NTT(ll *f,ll flag,ll n)
{
for(ll i=0;i<n;i++)
if(i<tr[i])swap(f[i],f[tr[i]]);
for(ll p=2;p<=n;p<<=1)
{
ll len=(p>>1);
ll yg=Quickly_Power(flag?G:invG,(mod-1)/p);
for(ll k=0;k<n;k+=p)
{
ll buf=1;
for(ll i=k;i<k+len;i++)
{
ll tmp=buf*f[i+len]%mod;
f[i+len]=ck(f[i]-tmp+mod);
f[i]=ck(f[i]+tmp),buf=buf*yg%mod;
}
}
}
if(!flag)
{
ll ny=Quickly_Power(n,mod-2);
for(ll i=0;i<n;i++)f
[i]=f[i]*ny%mod;
}
}
void Inv(ll *f,ll *g,ll m)
{
if(m==1)
{
g[0]=Quickly_Power(f[0],mod-2);
return;
}
//倍增
Inv(f,g,(m+1)>>1);
ll n=1;
while(n<(m<<1))n<<=1;
for(ll i=0;i<n;i++)
tr[i]=((tr[i>>1]>>1)|((i&1)?n>>1:0)),h[i]=f[i];
for(ll i=m;i<n;i++)h[i]=0;
NTT(h,1,n),NTT(g,1,n);
//求 h(x)(或g(x))
for(ll i=0;i<n;i++)
g[i]=(2-h[i]*g[i]%mod+mod)*g[i]%mod;
NTT(g,0,n);
for(ll i=m;i<n;i++)g[i]=0;
}
int main()
{
n=read();
for(ll i=0;i<n;i++)f[i]=read();
Inv(f,g,n);
for(ll i=0;i<n;i++)
printf("%lld ",g[i]);
return 0;
}
那么我们的多项式求逆能干什么呢?
我们知道,多项式乘除求逆都是基础工具。
我们需要通过这些工具,解决多项式对指开根、三角函数以及插值求值的问题。
多项式对数函数
我们通过求解多项式对数问题来初探乘法以及求逆的强大工具。

简单分析数据大小可以发现,朴素的做法根本无法支持非常大的数据。
因此我们还需多项式科技来解决。
拆开函数
这里需要一点微积分常识。
这里有个复合函数,比较难弄,因此可以将其求导然后积分来拆开。

代码实现
通过上述转化过程,我们不难想出程序如何写:

写代码的时候只有求导和求积这两个步骤是没有写过的。因此我放出代码。
求导:
void Differential(ll *f,ll *g,ll n)
{
for(ll i=1;i<n;i++)
g[i-1]=i*f[i]%mod;g[n-1]=0;
}
求积:
void Integral(ll *f,ll *g,ll n)
{
for(ll i=1;i<n;i++)
g[i]=f[i-1]*Quickly_Power(i,mod-2)%mod;
g[0]=0;
}
然后剩下的就非常好写了!
时间复杂度 O (n log n)。
#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
const ll mod=998244353;
#define G 3
#define MAXN 400040
#define invG 332748118
#define ck(x) ((x)>=mod?(x)-mod:(x))
inline ll read()
{
ll q=0;char ch=' ';
while(ch<'0'||ch>'9')ch=getchar();
while(ch>='0'&&ch<='9')q=(q<<3)+(q<<1)+ch-'0',ch=getchar();
return q;
}
ll n,tr[MAXN],f[MAXN],g[MAXN],w[MAXN];
ll ff[MAXN],_f[MAXN],_g[MAXN],h[MAXN];
inline ll Quickly_Power(ll a,ll b)
{
ll ans=1;
while(b)
{
if(b&1)ans=ans*a%mod;
a=a*a%mod,b>>=1;
}
return ans;
}
void NTT(ll *f,ll flag,ll n)
{
for(ll i=0;i<n;i++)
if(i<tr[i])swap(f[i],f[tr[i]]);
for(ll p=2;p<=n;p<<=1)
{
ll len=(p>>1);
ll wn=Quickly_Power(flag?G:invG,(mod-1)/p);
for(ll k=0;k<n;k+=p)
{
ll buf=1;
for(ll i=k;i<k+len;i++)
{
ll tmp=buf*f[i+len]%mod;
f[i+len]=ck(f[i]-tmp+mod);
f[i]=ck(f[i]+tmp),buf=buf*wn%mod;
}
}
}
}
void Inv(ll *f,ll *g,ll m)
{
if(m==1)
{
g[0]=Quickly_Power(f[0],mod-2);
return;
}
Inv(f,g,(m+1)>>1);
ll n=1;
while(n<(m<<1))n<<=1;
ll invn=Quickly_Power(n,mod-2);
for(ll i=0;i<n;i++)
tr[i]=((tr[i>>1]>>1)|((i&1)?n>>1:0)),w[i]=f[i];
for(ll i=m;i<n;i++)w[i]=0;
NTT(w,1,n),NTT(g,1,n);
for(ll i=0;i<n;i++)
g[i]=(2-w[i]*g[i]%mod+mod)*g[i]%mod;
NTT(g,0,n);
for(ll i=0;i<m;i++)
g[i]=invn*g[i]%mod;
for(ll i=m;i<n;i++)g[i]=0;
}
void Mul(ll *f,ll *g,ll *p,ll n,ll m)
{
m+=n,n=1;
while(n<m)n<<=1;
for(ll i=0;i<n;i++)
tr[i]=((tr[i>>1]>>1)|((i&1)?n>>1:0));
ll invn=Quickly_Power(n,mod-2);
NTT(f,1,n),NTT(g,1,n);
for(ll i=0;i<n;i++)
p[i]=f[i]*g[i]%mod;
NTT(p,0,n);
for(ll i=0;i<n;i++)
p[i]=p[i]*invn%mod;
}
void Differential(ll *f,ll *g,ll n)
{
for(ll i=1;i<n;i++)
g[i-1]=i*f[i]%mod;
g[n-1]=0;
}
void Integral(ll *f,ll *g,ll n)
{
for(ll i=1;i<n;i++)
g[i]=f[i-1]*Quickly_Power(i,mod-2)%mod;
g[0]=0;
}
void Ln(ll *f,ll *g,ll n)
{
Differential(f,ff,n);
Inv(f,_f,n);
Mul(ff,_f,_g,n,n);
Integral(_g,g,n);
}
int main()
{
n=read();
for(ll i=0;i<n;i++)f[i]=read();
Ln(f,g,n);
for(ll i=0;i<n;i++)
printf("%lld ",g[i]);
return 0;
}
总结与预告
本期教程就这么结束了,大家一定收获了很多!
我们学习了多项式求逆、多项式对数函数。
当然,多项式算法一定不止这么一点。
我们还会接触更多非常有用的多项式算法哦qwqq
下面一张图总结了一下算法、数据结构的分类。

总而言之,我会尽我所能向大家分享算法、数据结构,或者写工程的教程。
至于需要用到的一点点微积分常识,我也会另写教程。
我们下期再见!
欢迎常来我的博客Https://chtholly-oier.blog.luogu.org/