算法设计与分析第一次编程作业 1227. 快速傅里叶变换

题目描述

给定一个 n n n次多项式和一个 m m m次多项式,请求出 F ( x ) F(x) F(x) G ( x ) G(x) G(x)的乘积。

输入格式

第一行两个整数分别表示 F ( x ) F(x) F(x) G ( x ) G(x) G(x)的次数。
第二行 n + 1 n+1 n+1个整数,从低到高表示 F ( x ) F(x) F(x)的系数。
第三行 m + 1 m+1 m+1整数,从低到高表示 G ( x ) G(x) G(x)的系数。

输出格式

一行 n + m + 1 n+m+1 n+m+1个数字,从低到高表示 F ( x ) ⋅ G ( x ) F(x)\cdot G(x) F(x)G(x)的系数。

样例输入

    1 2
    1 2
    1 2 1

样例输出

    1 4 5 2

数据范围

对于30%的数据, n , m ≤ 3000 n, m \le 3000 n,m3000
对于100%的数据, n , m ≤ 200000 n, m \le 200000 n,m200000, 保证所有系数在0到9之间。

解决方案

#include <iostream>
#include <complex>
#include <cmath>
using namespace std;
#define cp complex<double>
#define PI acos(-1)
/***
这一版是可以通过的版本
声明:本题所有主要代码均为本人编写。
参考:
1、PI的定义,采用了https://blog.csdn.net/weijifen000/article/details/82691555?ops_request_misc=%257B%2522request%255Fid%2522%253A%2522164653173016780265417960%2522%252C%2522scm%2522%253A%252220140713.130102334..%2522%257D&request_id=164653173016780265417960&biz_id=0&utm_medium=distribute.pc_search_result.none-task-blog-2~all~baidu_landing_v2~default-2-82691555.es_vector_control_group&utm_term=C%2B%2BPI&spm=1018.2226.3001.4187
2、主要思路采用了老师上课讲授的算法,在实现上做了一些改动。
3、复数类的使用,参考了:https://horizonwd.blog.csdn.net/article/details/81478582
致谢:非常感谢老师给我指出浮点数比较的问题,也感谢助教老师的关心。
***/
void fft(cp *pq, cp *pq_omega, cp omega, int D,int inv)
{
    //注意浮点数的比较,原来老师提示的abs(omega.real()-1)<1e-8仍然有误差。
    //之前自己写的omega.real()==1非常错误,会导致死循环
    if (D == 1) { 
        for (int i=0;i<D;++i) pq_omega[i] = pq[i];
        return;
    }
    cp *pq_e = new cp[D/2];
    cp *pq_o = new cp[D/2];
    cp *pq_e_omega = new cp[D/2];
    cp *pq_o_omega = new cp[D/2];
    for (int i=0;i<D/2;++i) {pq_e[i] = pq[2*i];}
    for (int i=0;i<D/2;++i) {pq_o[i] = pq[2*i+1];}

    fft(pq_e,pq_e_omega,omega*omega,D/2,inv);
    fft(pq_o,pq_o_omega,omega*omega,D/2,inv);
    pq_omega[0] = pq_e_omega[0]+pq_o_omega[0];//0是特殊情况
    pq_omega[D/2] = pq_e_omega[0]-pq_o_omega[0];
    for (int t=1;t<D/2;++t) {
    //为了简化omega^t的计算,此处利用inv和w的性质(注意,此处的D是每次调用fft时候产生的除以2以后的D)
        pq_omega[t] = pq_e_omega[t]+cp(cos(2*t*PI/D),inv*sin(2*t*PI/D))*pq_o_omega[t];
        pq_omega[t+D/2] = pq_e_omega[t]-cp(cos(2*t*PI/D),inv*sin(2*t*PI/D))*pq_o_omega[t]; 
    }
    delete pq_e;
    delete pq_o;
    delete pq_e_omega;
    delete pq_o_omega;
}
void multiply(cp *p, cp *q, cp *r, int D, int p_order, int q_order)
{
    cp omega(cos(2*PI/D),sin(2*PI/D));
    cp *p_omega = new cp[D];
    cp *q_omega = new cp[D];
    cp *r_omega = new cp[D];
    fft(p,p_omega,omega,D,1);
    fft(q,q_omega,omega,D,1);
    for (int t=0;t<D;++t) r_omega[t] = p_omega[t]*q_omega[t];
    fft(r_omega,r,cp(cos(2*PI/D),-sin(2*PI/D)),D,-1);
    for (int i=0;i<D;++i) r[i] /= D;
    delete p_omega;
    delete q_omega;
    delete r_omega;
}
int main()
{
    int D_p=0, D_q=0,p_order,q_order;
    cin >> p_order >> q_order;
    while ((1<<D_p)<((p_order+1)*2)) ++D_p;
    while ((1<<D_q)<((q_order+1)*2)) ++D_q;
    int D = (D_p > D_q ? D_p : D_q);
    D = 1 << D;

    int *p = new int[D];
    int *q = new int[D];
    cp *cp_p = new cp[D];
    cp *cp_q = new cp[D];
    cp *r = new cp[D];
    for (int i=0;i<p_order+1;++i) {cin >> p[i];cp_p[i] = cp(p[i],0);}
    for (int i=p_order+1;i<D;++i) {p[i] = 0;cp_p[i] = cp(0,0);}
    for (int i=0;i<q_order+1;++i) {cin >> q[i];cp_q[i] = cp(q[i],0);}
    for (int i=q_order+1;i<D;++i) {q[i] = 0;cp_q[i] = cp(0,0);}
    multiply(cp_p,cp_q,r,D, p_order, q_order);
    for (int i=0;i<p_order+q_order+1;++i) cout <<  (int)(r[i].real()+0.5) << ' ';
    delete p;
    delete q;
    delete cp_p;
    delete cp_q;
    delete r;
    return 0;
}

总结

本题是根据课上知识实现FFT的一个版本,采用了动态数组,会有点笨拙,但是基本上还原了算法。需要注意的点是:1、精度问题,在fft函数第一行的截断上,要用整数才能够精确地实现;2、fft函数后半部分的数学实现比较复杂,不能搞混。

声明

本人对代码保有版权,但是题目来自校内OJ平台,本人没有题目的版权。如有侵权,请联系本人删除。

  • 2
    点赞
  • 0
    收藏
    觉得还不错? 一键收藏
  • 0
    评论

“相关推荐”对你有帮助么?

  • 非常没帮助
  • 没帮助
  • 一般
  • 有帮助
  • 非常有帮助
提交
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

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

抵扣说明:

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

余额充值