HDU1402(FFT入门)
题目链接:http://acm.hdu.edu.cn/status.php?user=Reykjavik11207&pid=1402&status=5
本题数据范围为5e4,常规方法O(n2)肯定是不行的。
FFT是离散傅里叶变换DFT的快速形式
对多项式f(x) = a0 + a1x + a2x2 + ··· +an-1xn-1,有两种表示法:
系数表达式 : (a0 , a1 , ··· , an-1)
由于n-1次多项式需要n个点来确定
所以可以用点值表达式 : ( (x0,f(x0)) , (x1,f(x1)) , ··· , (xn-1,f(xn-1)) ) 来表示
要获得点值表达式,首先选取n个x值获得对应f(x)的值
将f(x)分为奇偶两个部分f(x) = a0 + a2x2 + ··· + an-2xn-2 + a1x + a3x3 + ··· + an-1xn-1,
令f1(x) = a0 + a2x + ··· + an-2 x(n-2)/2,
f2(x) = a1 + a3x + ··· + an-1 x(n-1)/2,
则有f(x) = f1(x2) + xf2(x2)
f1(x)与f2(x)再分别分解,直至到常数ai为止
到了这里,复杂度并没有降低,反而由于x的整次幂未知还升高了,可以发现x = 1可以使式子更简单,因为1的多少次幂都是1,然后就是-1,但只有2个数远远不够,所以引入了复数。
是复平面单位圆上逆时针按k从小到大均匀分布的复根,间隔角度为2π/n,所以有:
= cos(2kπ/n) + i * sin(2kπ/n) ,计算复根的k次幂显然较实数更为方便,(但STL中三角函数也不是O(1),是多少我也不太懂,总之我把wnk函数放三重循环里就超时了)。
容易看出它的周期为n,即满足
同时还有以下性质: = cos(2kπ/n + π) + i * sin(2kπ/n + π) =
= cos(2*2k*π/n) + i * sin(2*2k*π/n) = cos(2kπ/(n/2)) + i * sin(2kπ/(n/2)) =
故f(wnk+n/2) = f1(wn/2k) - wnk f2(wn/2k)
以n = 4为例:
显然n必须为2的整数次幂
原来第i个元素的位置经过变换后的位置为i的二进制按长度翻转,
如上图中(0)10 = (00)2 翻转 -> (00)2 = (0)10,
(1)10 = (01)2 翻转 -> (10)2 = (2)10 ,
(2)10 = (10)2 翻转-> (01)2 = (1)10 ,
(3)10 = (11)2 翻转-> (11)2 = (3)10 ,
然后自底向上代入函数中,得到n个f(x)值
另一个数同理得到n个g(x)值
令F(x) = f(x)*g(x),则F(xi) = f(xi)*g(xi)
接下来就要用傅里叶逆变换IDFT求出F(x)的系数表达式
变换矩阵
的逆矩阵为
进而得到F(x)的系数表达式,即为结果
#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
const int N = 1 << 17;
const double pi = acos(-1.0);
const double eps = 1e-4;//这个...精度问题,本来用的1e-8就WA了
int n,len;
int rev(int x)//二进制翻转
{int res = 0;for(int i = 0; i < len; i++)res += ((x >> i) & 1) << (len - 1 - i);return res;
}
struct Complex
{double real,image;Complex(double r = 0,double i = 0){real = r;image = i;}Complex operator + (const Complex &t){return Complex(real + t.real, image + t.image);}Complex operator - (const Complex &t){return Complex(real - t.real, image - t.image);}Complex operator * (const Complex &t){return Complex(real * t.real - image * t.image, real * t.image + t.real * image);}
};
Complex wnk(double n,double k)
{return Complex(cos(2*pi*k/n), sin(2*pi*k/n));
}
void fft(Complex y[], int dft)
{Complex t1,t2;for(int i = 1; i < n; i <<= 1){Complex W(1,0), w = wnk(2*i,dft);for(int k = 0; k < i; k++){for(int j = k; j < n; j += i<<1){t1 = y[j] + W * y[j+i];t2 = y[j] - W * y[j+i];y[j] = t1;y[j+i] = t2;}W = W * w;}}if(dft == -1){for(int i = 0; i < n; i++)y[i].real /= n;}
}
Complex a1[N],a2[N],a[N];
int ans[N];
char stra[N>>1],strb[N>>1];
int main()
{while(~scanf("%s %s",stra,strb)){int lena = strlen(stra);int lenb = strlen(strb);len = log10(lena+lenb)/log10(2) + 1 + eps;n = 1 << len;for(int i = 0; i < lena; i++)a1[rev(i)] = Complex((double)(stra[lena-i-1] - '0'), 0);for(int i = lena; i < n; i++)a1[rev(i)] = Complex(0,0);for(int i = 0; i < lenb; i++)a2[rev(i)] = Complex((double)(strb[lenb-i-1] - '0'), 0);for(int i = lenb; i < n; i++)a2[rev(i)] = Complex(0,0);fft(a1,1);fft(a2,1);for(int i = 0; i < n; i++)a[rev(i)] = a1[i] * a2[i];fft(a,-1);int t = 0;for(int i = 0; i < n; i++){ans[n - 1 - i] = ((int)(a[i].real + eps) + t) % 10;t = ((int)(a[i].real + eps) + t) / 10;}bool flag = 0;for(int i = 0; i < n - 1; i++){if(ans[i])flag = 1;if(!flag)continue;printf("%d",ans[i]);}printf("%d\n",ans[n-1]);}return 0;
}
⎢⎢⎢(W0n)0(W1n)0⋮(Wn−1n)0(W0n)1(W1n)1⋮(Wn−1n)1⋯⋯⋱⋯(W0n)n−1(W1n)n−1⋮(Wn−1n)n−1⎤⎦⎥⎥⎥⎥⎥⎡⎣⎢⎢⎢⎢a0a1⋮an−1⎤⎦⎥⎥⎥⎥
⎢⎢⎢(W0n)0(W1n)0⋮(Wn−1n)0(W0n)1(W1n)1⋮(Wn−1n)1⋯⋯⋱⋯(W0n)n−1(W1n)n−1⋮(Wn−1n)n−1⎤⎦⎥⎥⎥⎥⎥⎡⎣⎢⎢⎢⎢a0a1⋮an−1⎤⎦⎥⎥⎥⎥
转载于:https://www.cnblogs.com/westwind1005/p/8886973.html
HDU1402(FFT入门)相关推荐
- FFT快速傅里叶变换 超详细的入门学习总结
FFT快速傅里叶变换 说明 本文创作的目的是为自己巩固该算法,加深印象并深入理解,同时也为FFT入门学者提供一份可鉴的学习总结. 原文链接:https://blog.csdn.net/qq_39565 ...
- 洛谷P3338:力(FFT)
传送门 解析 算是比较适合的FFT入门题了吧 一个重要的trick: 当函数无法表示成卷积时,可以把函数翻转过来 然后调一调就又是卷积了 一个重要的注意事项是FFT的lim一定是两多项式相乘结果多项式 ...
- maomao的每日动向
\(2019.02.04\) \(Nothing\) \(to\) \(do\). \(2019.02.05\) - 早上睡到\(12\)点 - 中午下午:吃饭串门拜年 - 晚上:吹爆<流浪地球 ...
- OpenGL入门学习笔记(一)——简单实现FFT海洋
一.前言 文章不赘述OpenGL的使用入门,使用入门请参考LearnOpenGL CN(https://learnopengl-cn.github.io/). 文章主要参考: [1][学习笔记]Uni ...
- simulink中的FFT 小白入门
simulink中的FFT 小白入门 前言 一.FFT是什么? 二.simulink中的FFT 1.FFT的实现方式 2.具体操作步骤 1.搭建模型 2.运行仿真 3.FFT分析 总结 前言 simu ...
- HDU1402 A * B Problem Plus(FFT)
http://acm.hdu.edu.cn/showproblem.php?pid=1402 初学FFT. http://www.cnblogs.com/WABoss/p/FFT_Note.html ...
- FFT:从入门到沉迷
终于学会了FFTFFTFFT并无法自拔 啥是FFT FFTFFTFFT指快速傅里叶变换 在OIOIOI中的应用是O(nlogn)O(nlogn)O(nlogn)计算函数卷积 人话:多项式乘法 多项式毒 ...
- 在quartus中使用FFT IP核最全教程(从入门到放弃)
一.准备工作 首先需要把需要的器材准备好,我使用的是quartus18.0,并且要使用IP核被破解的版本,不然无法使用其中的FFT和NCO,一定要注意,quartus对于版本非常敏感,一定要严格对应好 ...
- FFT从入门到使用(ACM/OI)
(首先膜一发高中生qwq亏我学了信号处理原理专业课了竟然是看高中生的博客才能看懂qwq...大概断断续续看了几十篇吧... 序言 可能是许久以来写的最认真的一篇了吧,加上今天感觉不太舒服,但是觉得要开 ...
最新文章
- 脑内世界模型:脑科学基础上的意识问题哲学解说
- python 菜鸟入门
- 中的stop_谈谈stop容器
- 机器学习面试——XGBoost,GBDT,RF(上)
- Python深入-Python的内存管理
- spring加载xml配置文件
- ExtAspNet v2.0.6发布 - AJAX性能提升
- python写透视挂_python opencv 透视变换
- 190226每日一句
- ActionScript Adobe Flash Builder Adobe Flash CC 学习笔记
- 【RK2206】3. 处理迪文屏事件
- esp32运行linux,ubuntu系统搭建ESP32 开发环境
- 苏宁回应股权质押给淘宝;日本政府用 AI 帮民众找对象;魅族回应 “暗中给手机植入木马” | EA周报...
- UniPro、Bugzilla和Teambition 缺陷管理工具优劣势对比
- 如果宁静是Oracle,那万茜,张雨绮,黄圣依 是什么?(附姐姐信息表)
- 7-20 打印九九口诀表
- 大白菜装机教程win10_如何用光盘快速重装系统(图文教程)
- 【本人秃顶程序员】程序员不要去这样的公司
- 用python搭建 百万答题 、自动百度搜索答案
- 用1、3、5、7 这4 个数字,能组成的互不相同且无重复数字的三位数有哪些?共有多少个?这些数的和为多少?
热门文章
- 安装wampserver及配置php,phpmyadmin遇到的问题及解决方法
- 牛腩44 整合登陆页 RequiredFieldValidator 和 ValidationSummary 以及 asp.net 自带的MD5 加密...
- 挨踢人生路--记我的10年18家工作经历 - 后记
- python 微服务框架 知乎_序: 我需要一个什么样的微服务框架
- 软件性能测试vu脚本录制,利用LR插件完成性能测试脚本
- main方法_错误: 在类 ZiFUChuan.Pyramid 中找不到 main 方法, 请将 main 方法定义为:
- python 非线性回归_机器学习入门之菜鸟之路——机器学习之非线性回归个人理解及python实现...
- 博途v13打开软件时显示连接不到服务器,TIA博途V13软件在打开程序过程中出现以下这种情况,怎么回事?...
- 查询两张表 然后把数据并在一起_工作表数据查询时,类似筛选功能LIKE和NOT LIKE的应用...
- golang 字符串排序_Golang操作数据库Redis