博客
关于我
强烈建议你试试无所不能的chatGPT,快点击我
【FFT&NTT 总结】
阅读量:4588 次
发布时间:2019-06-09

本文共 4173 字,大约阅读时间需要 13 分钟。

 

$FFT$总结

(因为还不会啊,,都没做过什么题,所以一边学一边打咯。。

 

1、主要是用来加速卷积形式的求和吧?

  $F*G(n)=F[i] × G[n-i]$

  平时是$O(n^2)$的,FFT可以$O(nlogn)$

2、相当于求两个多项式的乘积(你要求的函数是其系数)

  $A(x)=A0+A1*x+A2*x^2+A3*x^3+...+An−1*x^{n−1}$

  $B(x)=B0+B1*x+B2*x^2+B3*x^3+...+Bm−1*x^{m−1}$

3、具体步骤?

  系数表达->点值表达->相乘->点值表达->系数表达

4、点值表示法

把多项式$A(x)$代入若干个x的值得到若干个点$(x0,A(x0)),(x1,A(x1)),(x2,A(x2)),...,(xn−1,A(xn−1))$

 

  我们把从点值表达式转化为系数表达式的操作称为插值

 [点值表示法对应系数表示法是有唯一性的]

5、n次单位复根

  n次单位复根是满足w^n=1的复数w,有n个。他们均匀分布在以复平面的原点为圆心的单位圆上

  为$e^{\dfrac{2πki}{n}}$

  复数幂的定义$e^{ui}=cos(u)+isin(u)$

 

不如看

很详细的。

于是略过。

 

递归式模板

 

#include
#include
#include
#include
#define N 400010using namespace std;const double pi=acos(-1); struct P{ double x,y; P() {x=y=0;} P(double x,double y):x(x),y(y){}}a[N],b[N]; P operator + (P x,P y) {return P(x.x+y.x,x.y+y.y);}P operator - (P x,P y) {return P(x.x-y.x,x.y-y.y);}P operator * (P x,P y) {return P(x.x*y.x-x.y*y.y,x.x*y.y+x.y*y.x);} void fft(P *s,int n,int t){ if(n==1) return; P a0[n>>1],a1[n>>1]; for(int i=0;i<=n;i+=2) a0[i>>1]=s[i],a1[i>>1]=s[i+1]; fft(a0,n>>1,t);fft(a1,n>>1,t); P wn(cos(2*pi/n),t*sin(2*pi/n)),w(1,0); for(int i=0;i<(n>>1);i++,w=w*wn) s[i]=a0[i]+w*a1[i],s[i+(n>>1)]=a0[i]-w*a1[i]; //w^2=(w+(n>>1))^2 均匀分布在圆上面? //w[i^2,n]=w[i/2,n/2] 折半引理 //s[i]=a0’(i^2)+i*a1’(i^2)=a0(i)+i*a1(i) //s[i+n>>1]=a0’((i+n>>1)^2)+i*a1’((i+n>>1)^2)=a0’(i^2)-i*a1’(i^2) //因为i=-(i+n>>1) 折半引理} int main(){ int n,m,nn; scanf("%d%d",&n,&m); memset(a,0,sizeof(a));memset(b,0,sizeof(b)); for(int i=0;i<=n;i++) scanf("%lf",&a[i].x); for(int i=0;i<=m;i++) scanf("%lf",&b[i].x); nn=1;while (nn<=n+m) nn<<=1; fft(a,nn,1);fft(b,nn,1); for(int i=0;i<=nn;i++) a[i]=a[i]*b[i]; fft(a,nn,-1); for(int i=0;i<=n+m;i++) printf("%d ",(int)(a[i].x/nn+0.5)); return 0;}

 

  

迭代式模板

高效实现$FFT$

也就是使用迭代代替递归

比如n=8,R数组长这样:

$a0,a4,a2,a6,a1,a5,a3,a7$

也称作位逆序置换

 

#include
#include
#include
#include
#include
#include
#include
#define Maxn 262145#define pi acos(-1)using namespace std;struct P{ double x,y; P() {x=y=0;} P(double x,double y):x(x),y(y){} friend P operator + (P x,P y) {return P(x.x+y.x,x.y+y.y);} friend P operator - (P x,P y) {return P(x.x-y.x,x.y-y.y);} friend P operator * (P x,P y) {return P(x.x*y.x-x.y*y.y,x.x*y.y+x.y*y.x);}}a[Maxn],b[Maxn];int nn;int R[Maxn];void fft(P *a,int f){ for(int i=0;i
>1]>>1)|((i&1)<<(ll-1)); fft(a,1);fft(b,1); for(int i=0;i<=nn;i++) a[i]=a[i]*b[i]; fft(a,-1); for(int i=0;i<=n+m;i++) printf("%d ",(int)(a[i].x/nn+0.5)); return 0;}

  

 


 

 

$NTT$

 

 

 

简单复习一下原根

  1、模p下原根g,即g在模p下的阶为$\varphi(p)$

  2、g的幂构成模p下的缩系。【很有用!

  3、p有原根当且仅当 p=2,4,质数^a,2*质数^a

  4、求原根的方法:

    暴力枚举判断,先找出最小的一个。

    判断方法:对$\varphi(p)$质因数分解,假设为$p1^{r1}*p2^{r2}...*pn^{rn}$

         有恒有 $g^{\dfrac{\varphi(p)}{pi}}!=1 (Mod P)$成立,则g是p的原根,否则不是。

    假设最小原根是$g$,则当$gcd(d,\varphi(p))==1$时,$g^d$也是模p下的原根。

    即模p下原根个数为$\varphi(\varphi(p))$

 

对于$NTT$的题目,原根代替了n次单位复根,作用大致相同。只需把一些地方改成Mod之类的就可以。

 

 

具体看 【这篇写得太好了

 

 

 

 

当模数不符合要求?

 

 

具体实现方式(不完整代码,只放主要部分)

#include
#include
#include
#include
#include
using namespace std;#define Maxn 500010#define LL long longconst int Mod=998244353;const int G=3;int nn,R[Maxn],inv;void ntt(int *s,int f){ for(int i=0;i
>1]>>1)|((i&1)<<(ll-1)); inv=qpow(nn,Mod-2); ntt(a,1);ntt(b,1); for(int i=0;i<=nn;i++) a[i]=1LL*a[i]*b[i]%Mod; ntt(a,-1); /// return 0;}

  

好啦!!!

其实有很多地方是要证明的,但是我都不会 先记住吧。。


 

 

 

 

 

 

多项式求逆的代码实现

void get_inv(int *a,int *b,int len)  {    static int temp[Maxn];    if(len==1)      {          b[0]=qpow(a[0],Mod-2);          b[1]=0;          return;      }      get_inv(a,b,len>>1);      memcpy(temp,a,sizeof(int)*len);      memset(temp+len,0,sizeof(int)*len);      NTT(temp,len<<1,1),NTT(b,len<<1,1);      for(int i=0;i<(len<<1);i++)        b[i]=1LL*b[i]*(2-1LL*temp[i]*b[i]%Mod+Mod)%Mod;      NTT(b,len<<1,-1);    memset(b+len,0,sizeof(int)*len);  }

 

 

len为当前模数。不断二分。

 

 

2017-04-14 15:00:52

转载于:https://www.cnblogs.com/Konjakmoyu/p/6708101.html

你可能感兴趣的文章
NET设计规范(二) 命名规范
查看>>
VMware 9.0.1安装Mac OS X Mountain Lion 10.8.2
查看>>
SSL延迟
查看>>
android新手关于左右滑动的问题,布局把<android.support.v4.view.ViewPager/><ImageView/> 放在上面就不行了。...
查看>>
深入理解DIP、IoC、DI以及IoC容器
查看>>
赋值文件
查看>>
Vue 数组 字典 template v-for 的使用
查看>>
蓝牙模块选择经验谈
查看>>
java中==和equals
查看>>
CCActionPageTurn3D
查看>>
python random
查看>>
esp32-智能语音-cli(调试交互命令)
查看>>
netty与MQ使用心得
查看>>
apache反向代解决绝对路径可能出现的问题
查看>>
Oracle Metadata
查看>>
jquery 实现3d切割轮播图
查看>>
学习spring cloud 笔记
查看>>
字符串截取,SubString
查看>>
Android: 网络随时需要在3G和Wifi切换,网络程序需要注意
查看>>
ajax调用servlet
查看>>