FFT temp

FFT

在多项式乘法中的应用。请添加图片描述

  • 如果用传统的系数表示法来表示多项式(比如:f(x) = x^2 + 2x + 1,g(x) = 2x + 1,求f(x) * g(x)), 那么最暴力的做法是两个for循环,遍历多项式的每一项并相乘,这样做的时间复杂度是O(n^2)。

  • 前置知识1:多项式还有另一种表达方式:点值表示法。也就是长度为n的多项式需要n个不同点,就能确定唯一表示(举个例子:(1,2), (3,4) 两点确定一条直线)

    比如有个多项式f(x) = x^2 + 2x + 1,将x = 0, x = 1, x = 2分别代入,就可以得到三个点(0, 1), (1, 4), (2, 9),通过这三个点就可以确定多项式各项的系数。
    如果知道了两个多项式的点值表示,那么在O(n)的时间内将他们分别相乘,就可以得到两个多项式相乘得到多项式的点值表示,最后由点值表示就可以确定多项式各项的系数。
    在这里插入图片描述

  • 前置知识2:数学上,n次单位根是n次幂为1的复数,它们位于复平面的单位圆上。
    请添加图片描述
    每一个n次单位根可以通过欧拉公式算出来。(因为复数的相乘为模长相乘,幅角相加)

(前置知识有点多。。感觉放组会上讲不太合适。。有讲课而不是讨论的感觉

(看到某位大佬写的很好,就直接盗上来了https://www.cnblogs.com/flashhu/p/9323521.html
在这里插入图片描述
DFT把多项式由系数表达式转化为了点值表达式,然后就可以两个点值表达式相乘,在通过IDFT转化成系数表达式,就得到多项式相乘后每个项的系数。
在这里插入图片描述

#include <bits/stdc++.h>
using namespace std;
const int N = 4.2e6;
const double PI = acos(-1);
complex<double> f[N],g[N];
void FFT(complex<double> *a, int n, int op){//同上
    if(!n) return;
    complex<double> a0[n], a1[n];
    for(int i = 0; i < n; ++i)
        a0[i] = a[i<<1], a1[i] = a[i<<1|1];
    FFT(a0, n>>1, op); FFT(a1, n>>1, op);
    complex<double>W(cos(PI/n), sin(PI/n) * op), w(1, 0);
    for(int i = 0; i < n; ++i, w *= W)
        a[i] = a0[i] + w * a1[i], a[i+n] = a0[i] - w * a1[i];
}
int main(){
    ios::sync_with_stdio(false); cin.tie(0);

    int n,m;
    scanf("%d%d",&n,&m);
    for(int i = 0; i <= n; ++i) cin >> f[i];
    for(int i = 0; i <= m; ++i) cin >> g[i];
    for(m + =n,n = 1; n <= m; n <<= 1);//把长度补到2的幂,不必担心高次项的系数
    FFT(f, n>>1, 1);FFT(g, n>>1, 1);//DFT
    for(i = 0; i < n; ++i) f[i] *= g[i];//点值直接乘
    FFT(f, n>>1, -1);//IDFT
    for(i = 0; i <= m; ++i) printf("%.0f ",fabs(f[i].real()) / n);
    puts("");
    return 0;
}
解锁代码#include"stdio.h" #include"math.h" #include"stm32f10x.h" typedef struct complex COMPLEX; //创建了一个类型别名 COMPLEX,它指向 complex 结构体,便于后续代码的编写 void fft(COMPLEX *x); //FFT int period(COMPLEX *T,int k); //周期计算程序 struct complex //构造复数结构 //这定义了一个名为 complex 的结构体,包含两个浮点数成员 real 和 imag,分别代表复数的实部和虚部 { float real; float imag; }; //FFT子程序 void fft(COMPLEX *x) //声明了一个名为 fft 的函数,该函数接受一个指向 COMPLEX 类型的数组的指针作为参数 { double temp_re; double temp_im; unsigned int i,j,k,P; unsigned int LH,K,B; COMPLEX tmp; LH=512; j=LH; //FFT函数定义 //初始化局部变量,其中 LH 是循环中的初始步长,j 是初始索引,tmp 用于交换复数 for(i=1;i<1023;i++) //倒序最高位加1,逢2向次高位进位 { if(i<j) { tmp.real=x[j].real; tmp.imag=x[j].imag; x[j].real =x[i].real; x[j].imag =x[i].imag; x[i].real=tmp.real; x[i].imag=tmp.imag; } K=512; while((j<K)==0) { j=j-K; K=K/2; } j=j+K; } //对1024点输入序列进行位 //实现了输入序列的位逆序排列,是FFT算法的一部分 for(i=1;i<=10;i++) //第一重循环控制蝶形的级数 { B=pow(2,(i-1)); for(j=0;j<=(B-1);j++) //控制每级的蝶形和旋转因子 { P=pow(2,(10-i))*j; temp_re=cos(2*3.1415926*P/1024); temp_im=sin(2*3.1415926*P/1024); for(k=j;k<=1023;k=k+pow(2,i)) { COMPLEX FFT_temp1,FFT_temp2; FFT_temp1.real=x[k].real; FFT_temp2.real=x[k+B].real; FFT_temp1.imag=x[k].imag; FFT_temp2.imag=x[k+B].imag; x[k].real=FFT_temp1.real+FFT_temp2.real*temp_re+FFT_temp2.imag*temp_im; x[k].imag=FFT_temp1.imag+FFT_temp2.imag*temp_re-FFT_temp2.real*temp_im; x[k+B].real=FFT_temp1.real-FFT_temp2.real*temp_re- FFT_temp2.imag*temp_im; x[k+B].imag=FFT_temp1.imag-FFT_temp2.i
03-08
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

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

抵扣说明:

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

余额充值