FFT
原理
(30条消息) 十分简明易懂的FFT(快速傅里叶变换)_路人黑的纸巾的博客-CSDN博客_fft
a(x)
是一个n次多项式,b(x)
是一个m次多项式要求两个多项式a(x),b(x)
的乘积h(x)
\(h(x)=a(x)*b(x)两个函数的卷积,直接求复杂度O(n^2)\),
两个函数卷积的傅里叶变化=两个函数傅里叶变化的乘积
\(F(a(x)*b(x))=F(a(x)) \times F(b(x))\)
将两个函数先进行傅里叶变化\(O(nlogn)\),在相乘\(O(n)\),在进行傅里叶逆变化\(O(nlogn)\),最终得到\(h(x)\)
模板如下
n,m是多项式的最高次,多项式有n+1,m+1的系数,从低到高排序.
#define MAXN 2097153
#define PI 3.141592653589
complex<double> f[MAXN], g[MAXN];
int n, m, rev[MAXN];
int bit = 0, len = 1;
void fft(complex<double>* a, int inv);
int main()
{
scanf("%d%d", &n, &m);
int x;
for (int i = 0; i <= n; i++)
{
scanf("%d", &x);
f[i].real(x);
}
for (int i = 0; i <= m; i++)
{
scanf("%d", &x);
g[i].real(x);
}
while (len < n + m + 2)
bit++, len *= 2;
fft(f, 1);
fft(g, 1);
for (int i = 0; i < len; i++)
f[i] = f[i] * g[i];
fft(f, -1);
for (int i = 0; i <= n + m; i++)
printf("%d ", int(f[i].real() / len + 0.5));
return 0;
}
void fft(complex<double>* a, int inv)
{
for (int i = 0; i < len; i++)
{
rev[i] = (rev[i >> 1] >> 1) | ((i & 1) << (bit - 1));
if (i < rev[i])swap(a[i], a[rev[i]]);
}
for (int mid = 1; mid < len; mid *= 2)
{
complex<double> temp(cos(PI / mid), inv * sin(PI / mid));
for (int i = 0; i < len; i += mid * 2)
{
complex<double> omega(1, 0);
for (int j = 0; j < mid; j++, omega *= temp)
{
complex<double> x = a[i + j], y = omega * a[i + j + mid];
a[i + j] = x + y, a[i + j + mid] = x - y;
}
}
}
}
例题:
P1919 【模板】A*B Problem 升级版(FFT 快速傅里叶变换) - 洛谷 | 计算机科学教育新生态 (luogu.com.cn)