首页 > 其他分享 >【学习笔记】快速傅里叶变换

【学习笔记】快速傅里叶变换

时间:2023-04-24 19:46:52浏览次数:40  
标签:const 变换 多项式 单位根 笔记 int 点值 cp 傅里叶

怎么有人省选后才来学FFT啊
由于时间原因,本篇笔记仅为个人总结,真正想要学习FFT的请参看这篇博客。


前置知识

单位根性质:

  • $ w_n^{2k}= w_{n/2}^k $
  • $ w_n^a +w_n^b =w_n^{a+b} $

算法原理

可知 n+1 个点可以唯一确定一条 n 次多项式,于是可以用 n 个点的点对集合表示一条曲线。
那么如果有两个用点值表示的多项式,将点值对应相乘就可得到两个多项式乘出的多项式。如果两个多项式的最高项次数分别为 n 、m,那么确定最终多项式的点需要 n+m+1 个。

设多项式的点值表示的序列为 \(F(x)\),那么令奇偶下标分别为一组,奇数为 \(Fl(x)\),偶数为 \(Fr(x)\),可得到:
\(F(x) = Fl(x^2) + x\times Fr(x^2)\)
如果将 x 以单位根替换,则:
\(F(w_n^k) = Fl(w_{n/2}^k) + w_n^k \times Fr(w_{n/2}^k)\)
\(F(w_n^{k+n/2}) = Fl(w_{n/2}^k) - w_n^k \times Fr(w_{n/2}^k)\)
可发现每次的数据范围是缩小了一半的,于是我们可以在 \(n \log n\) 的时间内分治求出多项式点值表示。
上述过程被称为 \(DFT\) 。

求出多项式的点值表示,很多时候我们需要的是原多项式及其系数,下列通过点值还原多项式的做法被称作 \(IDFT\)。
设计算出的点值序列为 \(G(x)\) ,则
$ G(x)=\sum_{i=0}^{n-1} F(x) ( w_n^i )^i $

$ n\times F(x) = \sum_{i=0}^{n-1} G(x) ( w_n^{-i} )^i $

证明将 1 式代入 2 式即可。
这其实相当于再做一遍\(DFT\),不过将单位根变成负的并将最终答案除二。

快速记忆:

F(x)=Fl(x)+x*Fr(x)
正负号考虑单位根的正负性,IDFT要除n

算法实现

上述实现有两个问题:A.常数巨大 B.精度问题
对于 B ,有更优秀的 \(NTT\) 来解决,现在暂且不提。
对于 A ,各路神仙对其有不同的常数优化方式。
首先单位根的计算是多次的,我们可以将递归展开,从底层向上遍历的方式减少计算。
然后数组的拷贝也是一大问题,下文介绍蝴蝶变换可以将数组顺序直接转换为最底层的顺序。

原来的递归版(数组下标,先偶后奇,从0开始):
0 1 2 3 4 5 6 7 第1层
0 2 4 6|1 3 5 7 第2层
0 4|2 6|1 5|3 7 第3层
0|4|2|6|1|5|3|7 第4层

我们要求的就是第四层的顺序。
发现一个数最终的位置就是它的二进制翻转,如 (6,110) 变为 (6,011),这也是可以用递归来做的。
代码更容易理解:

for(int i=0;i<n;i++)tr[i]=(tr[i>>1]>>1)|((i&1)?n>>1:0);

还有其他的优化什么的,建议参考各路神仙的博客。

luogu模板题代码
#include<bits/stdc++.h>
using namespace std;
inline int rd(){
	int f=1,j=0;
	char w=getchar();
	while(!isdigit(w)){
		if(w=='-')f=-1;
		w=getchar();
	}
	while(isdigit(w)){
		j=j*10+w-'0';
		w=getchar();
	}
	return f*j;
}
const int N=2000010,M=2350000;
const double pai=acos(-1);
int n,m;
struct cp{
	cp (double xx=0,double yy=0){x=xx,y=yy;}
	double x,y;
	cp operator +(cp const &b)const{return cp(x+b.x,y+b.y);}
	cp operator -(cp const &b)const{return cp(x-b.x,y-b.y);}
	cp operator *(cp const &b)const{return cp(x*b.x-y*b.y,x*b.y+y*b.x);}
}f[N*2],p[N*2];
int tr[N*2];
void fft(cp *f,bool flg){
	for(int i=0;i<n;i++)if(i<tr[i])swap(f[i],f[tr[i]]);
	for(int p=2;p<=n;p<<=1){
		int len=p>>1;
		cp tg(cos(2*pai/p),sin(2*pai/p));
		if(!flg)tg.y*=-1;
		for(int k=0;k<n;k+=p){
			cp buf(1,0);
			for(int l=k;l<k+len;l++){
				cp tt=buf*f[len+l];
				f[len+l]=f[l]-tt;
				f[l]=f[l]+tt;
				buf=buf*tg;
			}
		}
	}
	return ;
}
signed main(){
	n=rd(),m=rd();
	for(int i=0;i<=n;i++)f[i].x=rd();
	for(int i=0;i<=m;i++)p[i].x=rd();
	for(m+=n,n=1;n<=m;n<<=1);
	for(int i=0;i<n;i++)tr[i]=(tr[i>>1]>>1)|((i&1)?n>>1:0);
//	cout<<"kkk\n";
	fft(f,1),fft(p,1);
	for(int i=0;i<n;i++)f[i]=f[i]*p[i];
	fft(f,0);
	for(int i=0;i<=m;i++)printf("%d ",(int)(f[i].x/n+0.49));
	return 0;
}

标签:const,变换,多项式,单位根,笔记,int,点值,cp,傅里叶
From: https://www.cnblogs.com/flywatre/p/17350655.html

相关文章

  • BSGS(大步小步算法)学习笔记
    解决高次同余问题。\(a^x\equivb(\modp)\),其中\(a\)与\(p\)同余。这个形式与欧拉定理类似。思想:meetinthemiddle(折半搜索)。具体的,令\(x=A\timest-B\),且\(x\)一定在\([0,\phi(p))\)的范围内。但是\(p\)是质数时复杂度还是会爆炸。将\(x=A\timest-B\)带入......
  • Python学习笔记--json序列化时间报错-改源码
    问题:转换时间报错执行代码为:importjsonfromdatetimeimportdate,datetimed={"time1":date.today(),"time2":datetime.today()}res=json.dumps(d)#报错  TypeError:ObjectoftypedateisnotJSONserializable方案1:手动转换str()方案2:继承类......
  • Vue笔记汇总
    Vue3快速上手1.Vue3简介2020年9月18日,Vue.js发布3.0版本,代号:OnePiece(海贼王)耗时2年多、2600+次提交、30+个RFC、600+次PR、99位贡献者github上的tags地址:https://github.com/vuejs/vue-next/releases/tag/v3.0.02.Vue3带来了什么1.性能的提升打包大小减少41%初次......
  • TypeScript 学习笔记 — 数组常见的类型转换操作记录(十四)
    获取长度lengthtypeLengthOfTuple<Textendsany[]>=T["length"];typeA=LengthOfTuple<["B","F","E"]>;//3typeB=LengthOfTuple<[]>;//0取第一项FirstItemtypeFirstItem<Textendsany[]>......
  • 老杜Vue实战教程完整版笔记(一)Vue程序初体验
    Vue作为国内使用率最高,最火爆的前端框架学习这门技术也越来越重要~动力节点老杜最新版Vue2+3教程已经上线!还是原来的配方,还是熟悉的味道学习地址:https://www.bilibili.com/video/BV17h41137i4/1Vue程序初体验我们可以先不去了解Vue框架的发展历史、Vue框架有什么特点、Vu......
  • JAVA学习笔记随记1(类与对象)
    首先说明,这是为了学习java而做的笔记,所以记起来可能杂乱无章,无所谓了,刚开始学习都是这样的。。。首先小结下String的知识点String可以直接声明并赋初值并可以修改,例如:Stringabc="a";abc="b";其次字符串之间的连接用'+',只要出现字符串和其他数据类型之间用'+'连接,那么该......
  • VBA学习笔记
    2023-04-24(1)OptionExplicit:在模块最开始加这句代码,如果程序中有未声明的变量,程序不会运行,且计算机会自动提醒你声明变量。在VBE编辑界面,通过工具--》选项--》编辑器--》勾选“要求变量声明”,则每个模块都会在第一句自动写下“OptionExplicit”(2)Static关键字:声明变量为......
  • 前端进化笔记-JavaScript(一)
    简介:实现:三部分ecmascript语言核心DOM文档对象模型BOM浏览器对象模型DOM:用于html的应用程序接口(API),把整个页面映射成一个多层节点结构。例如:<html> <head> <title>samplepage</title> </head> <body> <p>helloworld</p> </body></html>DOM......
  • SpringMVC 框架总结笔记
    第一章初识SpringMVC1.1SpringMVC概述SpringMVC是Spring子框架SpringMVC是Spring为【展现层|表示层|表述层|控制层】提供的基于MVC设计理念的优秀的Web框架,是目前最主流的MVC框架。SpringMVC是非侵入式:可以使用注解让普通java对象,作为请求处理器【Controller】......
  • 学习笔记10
    第21章存储秘密21.1磁盘存储秘密的一个很直接的办法是把秘密存储在计算机的硬盘上或其他永久存储介质上,这是可行的,但是任何使用此电脑的人都能使用该密钥。一个更好的解决方案是让Alice把密钥存储在她的PDA或智能手机上。这些设备很少会借给别人使用,而且无论去哪里都会随......