感觉上次学习DFT简直是乱来了。不知道误导了多少人，这里深感抱歉。
这次我再看了看《算法导论》，觉得收获很大，终于粗略的知道DFT的原理了！
如何将两个多项式相乘
对于一个n次多项式，(f(x)=sum_{i=0}^n a_ix^i)，根据书上说的是有一个叫插值的东西。通俗的说，就是我们随意选取(n)个不同的(x)，不妨设为(x_0,x_1,...,x_n)，通过这个我们可以算出(y_i=f(x_i))，这个过程就是DFT。
我们现在需要将两个多项式相乘，所以这两个多项式都算出它们的(y_i)，分别记为(y_i^{[0]})和(y_i^{[1]})，然后我们得到(y_i=y_i^{[0]}y_i^{[1]})，然后将(y)进行IDFT就是答案了。根据插值多项式的唯一性，可以证明方法的正确性。
这里要注意的是，如果原来的两个多项式的项数分别为(n,m)，则我们要取(2^s(sin Z^+,2^s>n+m))个不同的(x)。
如何采用分治的方法计算DFT
如果我们要计算$$f(x)=sum_{i=0}{n-1}a_ixi$$
这里的(n=2^s(sin Z^+))。
设$$g(x)=sum_{i=0}{n/2-1}a_{2i}xi$$
那么(f(x)=g(x^2)+xh(x^2))。既然要分治，我们需要把问题的规模减少一半，但是(x^2)仍有(n)个取值。
为了减少(x^2)的取值，(x)的取值中我们可以有(n/2)对相反数，这样它们的平方是相同的。我们在确定这一轮(x)的取值后，下一轮(x)的取值就确定了。所以，如果我们简单的取(x=1,-1,2,-2,3,-3,...)，下一轮的(x=1,4,9,...)，就没有相反数了。
为了解决这个问题，单位复数根来了。
设(omega_n^k=e^{2pi ki /n})，(i)是虚数单位。
如果这一轮我们的(x=omega_n^k(k=0,1,2,...,n-1))，那么(omega_n^k)与(omega_n^{k+n/2}(k=0,1,2,...,n/2-1))将会是相反数。
而且，下一轮的(x=omega_{n/2}^{k}(k=0,1,2,...,n/2-1))将有(n/4)对相反数，问题规模恰好缩小为原来的一半。
如此，便可以得到以下代码：
fft(a)
n = a.length
if n == 1
return a
wn = cos(2 * PI / n) + sin(2 * PI / n)
w = 1
a0 = (a[0],a[2],a[4],...,a[n-2])
a1 = (a[1],a[3],a[5],...,a[n-1])
y0 = fft(a0)
y1 = fft(a1)
for k = 0 to n / 2 - 1
y[k] = y0[k] + w * y1[k]
y[k + n / 2] = y0[k] - w * y1[k]
w = w * wn
return y
那么IDFT？
我们上面的计算过程即是：
不会画矩阵，只好用表格替代了。
然后我们只要求出那个(n imes n)的矩阵的逆矩阵就可以了。
设原来的矩阵(V_n)的((k,j))处的元素为(omega_n^{kj})，可以证明，(V_n^{-1})的((j,k))处的元素为(omega_n^{-kj}/n)。
证明：考虑((j,j'))处的元素，$$sum_{k=0}{n-1}(omega_n{-kj}/n)(omega_n{kj'})=sum_{i=0}{n-1}omega_n^{k(j'-j)}/n$$
如果(j=j')，那么结果显然是(1)，否则，是(0)，因为它们都有自己的相反数，可以消掉。
这样就证明了它们的乘积是一个单位矩阵。
可以发现，我们只需要把原来的(omega_n^{kj})变为(omega_n^{-kj})就可以了，记得最后要除以(n)。这样的计算过程和原来的分治方法是一样的。
蝴蝶变换
其实这是一个对于OI来说不可忽视的常数优化。
我们来跟踪一下输入fft第(i)位的元素到了哪个位置，可以发现，如果(i=(11010))，那么最后到达的位置就是((01011))(这里括号里的是二进制)，就是把它的二进制翻转。于是就有下面的代码：
一开始要把元素移动到对应的位置。
void fft(cpd* A, bool rev) {
for (int i = 0; i < n; i ++)
if (bitRev[i] < i) swap(A[i], A[bitRev[i]]);
for (int s = 1; s <= MAXLOG; s ++) {
int m = 1 << s;
cpd unit(cos(PI * 2 / m), sin(PI * 2 / m) * (rev ? -1 : 1));
for (int k = 0; k + m <= n; k += m) {
cpd w(1, 0);
for (int j = 0; j < m >> 1; j ++) {
cpd t = w * A[k + j + (m >> 1)];
cpd u = A[k + j];
A[k + j] = u + t;
A[k + j + (m >> 1)] = u - t;
w *= unit;
}
}
}
if (rev)
for (int i = 0; i < n; i ++) A[i] /= n;
}
模下的
一般模数为(b2^a+1)而且是个质数，例子(7×17×2^{23}+1=998244353)。
同样的，我们要找到(n)个不同的(x)，也就是说，我们要找到一个(x)，使得(x^n=1)且(x^{n/2}=-1)，这里的等号都是模意义下的。
我们可以采用下面的方法找到这样的(x)：
因为(p)是质数，所以(a^{p-1}=1)（(a)是任意整数），不妨令(x=y^{(p-1)/n})，由于(n)是(2)的幂，所以这里一定有(n|p-1)。
接着我们就满足了(x^n=1)这一个条件，另外一个的话，我们可以枚举(y)，然后就可以找到了。
对于(998244353)，(y)可以取(3)。
记录下一个求素数(p)原根的方法：
由于原根比较小，我们可以从(2)开始枚举，现在的问题是如何判断一个数(a)是不是(p)的原根。
由于(a^{p-1}=1)，所以我们枚举(p-1)的素因子(f)，如果存在这样的(f)使得(a^{p-1 over f}=1)那么(a)就不是原根。