数学 (OI)

数学 (OI)

前言

说实话这首歌有点反。

Jay of Love

Another sunrise, another sunset
Soon it'll all be yesterday
Another good day, another bad day,
What did you do today?
Why do we choose to chase what we'll lose?
What you want isn't what you have.
What you have may not be yours, to keep.
If I could find love, at a stop, in a park with open arms,
I would save all my love, in a jar,
made of sparks, sealed in my beating heart,
Could it be yours to keep, the Jar of Love.
Another left turn, another head turns
Could he be someone I deserve?
Another right turn, another lesson learned
Never leave an open flame to burn
Why do we choose to chase what we'll lose?
What you want isn't what you have.
What you have may not be yours, to keep.
If I could find love,
at a stop, in a park with open arms,
I would save all my love,
in a jar,made of sparks,
sealed in my beating heart,
Could it be yours to keep, the Jar of Love.
Could you be my love
Could you be my love
Could you be my love
Could you be my love
Could you be her love
Could you be his love
Could you be my love
Could I be you love
If I could find love,
at a stop, in a park with open arms,
I would save all my love,
in a jar,made of sparks,
sealed in my beating heart,
Could it be yours to keep
If I could find love,
at a stop, in a park with open arms,
I would save all my love,
in a jar,made of sparks,
sealed in my beating heart,
Could it be yours to keep
If I could find love,
at a stop, in a park with open arms,
I would save all my love,
in a jar,made of sparks,
sealed in my beating heart,
Could it be yours to keep
the Jar of Love.
Could it be yours to keep
the Jar of Love.
Could it be yours to keep
the Jar of Love.

可能是最后一篇学习笔记。

FFT 学习笔记

傅里叶

很久之前就学了，现在好像忘了，连着 ntt,fwt 重新学一下。顺便打个学习笔记。

把 \(n\) 当做 \(2^k\)

用复数做多项式的根，具体的，\(x = w^1_n\) 一次单位根, 即 \((\cos \frac{2 \pi}{n}, \sin \frac{2 \pi}{n})\)。

理由的话可以证一下:

设 \(A(x) = \sum\limits_{i = 0}^{n - 1}a_ix^i\) 多项式，\(x\) 做 \(0 \dots n - 1\) 次单位根带出来的离散傅里叶变换为 \((y_0, \dots, y_{n - 1})\), 接着设 \(B(x) = \sum\limits_{i = 0}^{n - 1}y_ix^i\), 将 \(x^{-1}\) 带入作单位根得出新的离散傅里叶变换 \((z_0, \dots, z_{n - 1})\)。

那按题意模拟，

\[\begin{aligned}z_k &= \sum_{i = 0}^{n - 1}y_i \left(\omega_{n}^{k}\right)^{-i} \\&= \sum_{i = 0}^{n - 1}\left(\sum_{j = 0}^{n - 1}a_j \left(\omega_{n}^{j}\right)^j\right) \left(\omega_{n}^{k}\right)^{-i} \\&= \sum_{j = 0}^{n - 1} a_j\left(\sum_{i = 0}^{n - 1}(\omega_{n}^{j - k})^i\right) \\&= na_k \end{aligned}\]

所以这个玩意儿就是做两遍傅里叶变换的事儿。

快速

然后就是优化复杂度的问题了。傅里叶叔叔作为一个连电脑都没见过的人，他是闲得蛋疼才会优化复杂度，于是他蛋完好。

后人用分治优化了 fft 复杂度。

具体原理如下：

设 \(A(x) = a_0x^0 + \dotsb + a_{n - 1}x^{n - 1}\)， 按奇偶分个类：

\[A(x) = (a_0 + a_2x^2 + \dotsb + a_{n - 2}x^{n - 2}) + (a_1x + \dotsb + a_{n - 1}x^{n - 1}) \]

整俩多项式：

\[A_1(x) = a_0 + a_2x + \dotsb + a_{n - 2}x^{\frac{n}{2} - 1} \]

\[A_2(x) = a_1 + a_3x + \dotsb + a_{n - 1}x^{\frac{n}{2} - 1} \]

显显显 \(A(x) = A_1(x^2) + xA_2(x^2)\)。

假设 \(k < \frac{n}{2}\)，将 \(x = \omega_{n}^{k}\) 带入，

\[\begin{aligned} A(\omega_{n}^k) &= A_1(\omega_n^{2k}) + \omega_n^kA_2(\omega_n^{2k}) \\ &= A_1(\omega_{\frac{n}{2}}^{k}) + \omega_n^kA_2(\omega_{\frac{n}{2}}^{k}) \end{aligned} \]

对于 \(A(\omega_{n}^{k + \frac{n}{2}})\)：

\[\begin{aligned} A(\omega_{n}^{k + \frac{n}{2}}) &= A_1(\omega_n^{2k + n}) + \omega_n^{k + \frac{n}{2}}A_2(\omega_n^{2k + n}) \\ &= A_1(\omega_{\frac{n}{2}}^{k}) - \omega_n^kA_2(\omega_{\frac{n}{2}}^{k}) \end{aligned}\]

所以只要对于 \(A_1A_2\) 在 \((\omega_{\frac{n}{2}}^{0} \dotsb \omega_{\frac{n}{2}}^{\frac{n}{2} - 1})\) 就可以求出 \(A\) 在 \((\omega_{n}^{0} \dotsb \omega_{n}^{n - 1})\) 离散傅里叶变换。

\(n - 1\) 直接 return。

实现

其实还有递归版，就是常数忒他妈大，写了也没啥意义，就不写了，这里是倍增法，比较快。

然后就是各位 大牛 们纷纷拽高级名词，什么蝴蝶变换，根本看不懂。

然后看代码就是用个变量数组存起来，减少常数。

\(QwQ\)， 有点唐。

考虑一下非递归实现就是要从下往上走，那么我们就需要最终位置。

简单手摸一下得出， \(x\) 的最后位置是其二进制翻转的值。

一点一点向上还原即可。

CODE

#include <bits/stdc++.h>
typedef long long ll; 
using namespace std; 
const int N = (1 << 21) + 100; 
class _complex {
private:
double x, y; 
public:
_complex(double X = 0, double Y = 0): x(X), y(Y) {}
inline double Real() {return x; }
inline double Imag() {return y; }
inline _complex operator + (_complex a) {return _complex(x + a.x, y + a.y); }
inline _complex operator - (_complex a) {return _complex(x - a.x, y - a.y); }
inline _complex operator * (_complex a) {return _complex(x * a.x - y * a.y, a.x * y + x * a.y); }
} a[N], b[N]; 
int n, m, num = 1; 
const double Pi = acos(-1.0); 
int l, r[N]; void FFT(_complex *A, int type) {for (int i = 0; i < num; ++ i) {if (i < r[i]) swap(A[i], A[r[i]]); }for (int mid = 1; mid < num; mid <<= 1) {_complex Wn = _complex(cos(Pi / mid), sin(Pi / mid) * type); for (int R = mid << 1, j = 0; j < num; j += R) {_complex w = _complex(1, 0); for (int k = 0; k < mid; ++ k, w = w * Wn) {_complex x = A[j + k], y = w * A[j + k + mid]; A[j + k] = x + y, A[j + mid + k] = x - y; }}}
}signed main() {freopen("1.in", "r", stdin); freopen("1.out", "w", stdout); ios::sync_with_stdio(0), cin.tie(0), cout.tie(0); cin >> n >> m; for (int i = 0, x; i <= n; ++ i) cin >> x, a[i] = _complex(x, 0); for (int i = 0, x; i <= m; ++ i) cin >> x, b[i] = _complex(x, 0); while (num <= n + m) num <<= 1, ++ l; for (int i = 0; i < num; ++ i) r[i] = (r[i >> 1] >> 1) | ((i & 1) << (l - 1)); FFT(a, 1), FFT(b, 1); for (int i = 0; i < num; ++ i) a[i] = a[i] * b[i]; FFT(a, -1); for (int i = 0; i <= n + m; ++ i) cout << (int)(a[i].Real() / num + 0.5) << ' '; 
}