快速傅里叶变换不能三言两语能解释清楚,自己看了一些资料,仍不敢说完全掌握了。
快速傅里叶变换
(FFT)的作用及解释:
http://blog.jobbole.com/70549/
编程实现:
大神们的博客:
http://blog.miskcoo.com/2015/04/polynomial-multiplication-and-fast-fourier-transform#i-15
卷积是两个变量在某范围内相乘后求和的结果。如果卷积的变量是序列x(n)和h(n),则卷积的结果
,其中星号*表示卷积。
多项式乘法就可被FFT高效解决。
快速傅里叶变换提供了一种新的视觉,从不同的角度看待问题,让复杂问题简单化。把数字统统放到复数域内,由单位复根W的三大特性,对称性,可约性,周期性实现普通离散变换效率上的飞跃。
复根:
对称性奇偶分类,分治处理,此外
周期性:
可约性:
帮助减少了运算量。
快速傅里叶变换解决问题的过程可以说就是先将函数表达变成点值关系,进行运算,再转成系数表示。
A * B Problem Plus
http://acm.hdu.edu.cn/showproblem.php?pid=1402
Calculate A * B. 数据很大
#include <iostream>
#include <cstdio>
#include <cstring>
#include <cmath>
#include <algorithm>
using namespace std;
const int N=5e5+10;
const double PI=acos(-1.0);
struct comx{ // 复数
double r,i;
comx(double r=0.0,double i=0.0){
this->r=r;
this->i=i;
}
comx operator +(const comx x){
return comx(r+x.r,i+x.i);
}
comx operator -(const comx x){
return comx(r-x.r,i-x.i);
}
comx operator *(const comx x){
return comx(r*x.r-i*x.i,r*x.i+i*x.r);
}
};
// 倒位序·雷德算法 100 --> 001
void Rader(comx F[],int len){
int j=len/2;
for(int i=1;i<len-1;i++){
if(i<j) swap(F[i],F[j]);
int k=len>>1;
while(j>=k){
j-=k;
k>>=1;
}
if(j<k) j+=k;
}
}
// FFT 快速傅里叶变换
void FFT(comx F[],int len,int f){
Rader(F,len); // 倒位序 奇偶分类
for(int i=2;i<=len;i<<=1){ // 分治
comx I(cos(-f*2*PI/i),sin(-f*2*PI/i));
for(int j=0;j<len;j+=i){
comx w(1,0); // 旋转因子
for(int k=j;k<j+i/2;k++){
comx u=F[k];
comx t=w*F[k+i/2];
F[k]=u+t;
F[k+i/2]=u-t; //蝴蝶
w=w*I;
}
}
}
if(f==-1) {
for(int i=0;i<len;i++) F[i].r/=len;
}
}
// 卷积
void Conv(comx a[],comx b[],int len){
FFT(a,len,1); // DFT 变成点值
FFT(b,len,1);
for(int i=0;i<len;i++) a[i]=a[i]*b[i];
FFT(a,len,-1); // IDFT 变回关系式
}
char str1[N],str2[N];
comx a[N],b[N];
int ans[N],len;
int main()
{
while(~scanf("%s%s",str1,str2)){
int len1=strlen(str1),len2=strlen(str2);
len=1;
while(len<2*len1 || len<2*len2) len<<=1;
int i;
for(i=0;i<len1;i++){ //大数存入两个复数数组中
a[i].r=str1[len1-1-i]-'0';
a[i].i=0;
}
while(i<len){
a[i].r=a[i].i=0;
i++;
}
for(i=0;i<len2;i++){
b[i].r=str2[len2-1-i]-'0';
b[i].i=0;
}
while(i<len){
b[i].r=b[i].i=0;
i++;
}
Conv(a,b,len);
for(i=0;i<len;i++) ans[i]=a[i].r+0.5; //四舍五入
for(i=0;i<len;i++){
ans[i+1]+=ans[i]/10; //高精度处理
ans[i]%=10;
}
int length=0;
for(i=len-1;i>=0;i--){
if(ans[i]) {
length=i;
break;
}
}
for(i=length;i>=0;i--) printf("%d",ans[i]);
printf("\n");
}
return 0;
}
拓展:
将卷积部分改成单纯的求和,即成为高精度求和
求和 (前提是两个大数均是非负数):
void Conv(comx a[],comx b[],int len){
FFT(a,len,1); // DFT 变成点值
FFT(b,len,1);
for(int i=0;i<len;i++) a[i]=a[i]+b[i];
FFT(a,len,-1); // IDFT 变回关系式
}
但是不能拓展成为高精度减法:
void Conv(comx a[],comx b[],int len){
FFT(a,len,1); // DFT 变成点值
FFT(b,len,1);
for(int i=0;i<len;i++) a[i]=a[i]-b[i];
FFT(a,len,-1); // IDFT 变回关系式
}
POJ 2389 Bull Math
http://poj.org/problem?id=2389
普通的大数乘法。我就是试试FFT。。。
和上题差不多,数组大小改一下就能过
1028 大数乘法 V2
http://www.51nod.com/onlineJudge/questionCode.html#!problemId=1028
也是大数乘法
FFT用于多项式乘法,当数据很大时,其优势体现的更加明显,普通的计算过程为O(n^2),但是FFT可以做到O(nlogn).