注:考虑到对高中生友好的需要,本文尽量减少了对于高数和线代等前置知识的依赖。如果你是正在学习FFT的oi选手,可以放心食用。
我当年在学习FFT的时候,大部分文章都会在讲完DFT的过程之后说“对于IDFT,我们只需要把单位根全部取个共轭,最后把结果除以
n
n
n 就行了”之类的话,当然它们会对此做出解释:DFT的过程本质上就是在做一个线性变换,我们对着线性变换的矩阵求个逆折腾一下就能变换回去了。
设原始的多项式是
f
(
x
)
=
a
0
+
a
1
x
+
a
2
x
2
+
.
.
.
+
a
n
−
1
x
n
−
1
f(x) = a_0 + a_1x + a_2x^2 + ... + a_{n-1}x^{n-1}
f(x)=a0+a1x+a2x2+...+an−1xn−1 ,我们要对
k
=
0
,
1
,
.
.
.
,
n
−
1
k=0,1,...,n-1
k=0,1,...,n−1 求
b
k
=
f
(
ω
k
)
b_k=f(\omega^k)
bk=f(ωk) (其中
ω
=
e
2
π
i
n
\omega=e^{\frac{2\pi i}{n}}
ω=en2πi )。这一过程可以用矩阵表示为:
(
1
1
1
⋯
1
1
ω
ω
2
⋯
ω
n
−
1
1
ω
2
ω
4
⋯
ω
2
(
n
−
1
)
⋮
⋮
⋮
⋱
⋮
1
ω
n
−
1
ω
2
(
n
−
1
)
⋯
ω
(
n
−
1
)
2
)
(
a
0
a
1
a
2
⋮
a
n
−
1
)
=
(
b
0
b
1
b
2
⋮
b
n
−
1
)
\left( \begin{matrix} 1 & 1 & 1 & \cdots & 1\\ 1 & \omega & \omega^2 & \cdots & \omega^{n-1}\\ 1 & \omega^2 & \omega^4 & \cdots & \omega^{2(n-1)}\\ \vdots & \vdots & \vdots & \ddots & \vdots\\ 1 & \omega^{n-1} & \omega^{2(n-1)} & \cdots & \omega^{(n-1)^2}\\ \end{matrix} \right) \left( \begin{matrix} a_0\\ a_1\\ a_2\\ \vdots\\ a_{n-1} \end{matrix} \right)= \left( \begin{matrix} b_0\\ b_1\\ b_2\\ \vdots\\ b_{n-1} \end{matrix} \right)
⎝⎜⎜⎜⎜⎜⎛111⋮11ωω2⋮ωn−11ω2ω4⋮ω2(n−1)⋯⋯⋯⋱⋯1ωn−1ω2(n−1)⋮ω(n−1)2⎠⎟⎟⎟⎟⎟⎞⎝⎜⎜⎜⎜⎜⎛a0a1a2⋮an−1⎠⎟⎟⎟⎟⎟⎞=⎝⎜⎜⎜⎜⎜⎛b0b1b2⋮bn−1⎠⎟⎟⎟⎟⎟⎞
为方便起见,我们设
V
=
(
1
1
1
⋯
1
1
ω
ω
2
⋯
ω
n
−
1
1
ω
2
ω
4
⋯
ω
2
(
n
−
1
)
⋮
⋮
⋮
⋱
⋮
1
ω
n
−
1
ω
2
(
n
−
1
)
⋯
ω
(
n
−
1
)
2
)
,
α
=
(
a
0
a
1
a
2
⋮
a
n
−
1
)
,
β
=
(
b
0
b
1
b
2
⋮
b
n
−
1
)
V=\left( \begin{matrix} 1 & 1 & 1 & \cdots & 1\\ 1 & \omega & \omega^2 & \cdots & \omega^{n-1}\\ 1 & \omega^2 & \omega^4 & \cdots & \omega^{2(n-1)}\\ \vdots & \vdots & \vdots & \ddots & \vdots\\ 1 & \omega^{n-1} & \omega^{2(n-1)} & \cdots & \omega^{(n-1)^2}\\ \end{matrix} \right),\ \alpha= \left( \begin{matrix} a_0\\ a_1\\ a_2\\ \vdots\\ a_{n-1} \end{matrix} \right),\ \beta= \left( \begin{matrix} b_0\\ b_1\\ b_2\\ \vdots\\ b_{n-1} \end{matrix} \right)
V=⎝⎜⎜⎜⎜⎜⎛111⋮11ωω2⋮ωn−11ω2ω4⋮ω2(n−1)⋯⋯⋯⋱⋯1ωn−1ω2(n−1)⋮ω(n−1)2⎠⎟⎟⎟⎟⎟⎞, α=⎝⎜⎜⎜⎜⎜⎛a0a1a2⋮an−1⎠⎟⎟⎟⎟⎟⎞, β=⎝⎜⎜⎜⎜⎜⎛b0b1b2⋮bn−1⎠⎟⎟⎟⎟⎟⎞
其中
V
V
V 是
n
n
n 阶复矩阵,
α
,
β
\alpha,\beta
α,β是两个
n
n
n 维复向量,那么原式可以写作
V
α
=
β
V\alpha=\beta
Vα=β 。
DFT的任务是给定向量
α
\alpha
α 求向量
β
\beta
β ,IDFT的任务刚好相反,是给定向量
β
\beta
β 求向量
α
\alpha
α 。
等式两边同时左乘
V
−
1
V^{-1}
V−1 ,我们有
α
=
V
−
1
β
\alpha=V^{-1}\beta
α=V−1β,所以看起来我们只需要求出矩阵
V
V
V 的逆矩阵
V
−
1
V^{-1}
V−1 就可以了。这个时候,大部分文章都会直接告诉你:
V
−
1
=
1
n
(
1
1
1
⋯
1
1
ω
−
1
ω
−
2
⋯
ω
−
(
n
−
1
)
1
ω
−
2
ω
−
4
⋯
ω
−
2
(
n
−
1
)
⋮
⋮
⋮
⋱
⋮
1
ω
−
(
n
−
1
)
ω
−
2
(
n
−
1
)
⋯
ω
−
(
n
−
1
)
2
)
V^{-1}=\frac{1}{n}\left( \begin{matrix} 1 & 1 & 1 & \cdots & 1\\ 1 & \omega^{-1} & \omega^{-2} & \cdots & \omega^{-(n-1)}\\ 1 & \omega^{-2} & \omega^{-4} & \cdots & \omega^{-2(n-1)}\\ \vdots & \vdots & \vdots & \ddots & \vdots\\ 1 & \omega^{-(n-1)} & \omega^{-2(n-1)} & \cdots & \omega^{-(n-1)^2}\\ \end{matrix} \right)
V−1=n1⎝⎜⎜⎜⎜⎜⎛111⋮11ω−1ω−2⋮ω−(n−1)1ω−2ω−4⋮ω−2(n−1)⋯⋯⋯⋱⋯1ω−(n−1)ω−2(n−1)⋮ω−(n−1)2⎠⎟⎟⎟⎟⎟⎞
并且告诉你这一事实的正确性“不难验证”(也确实不难)。然后我们就可以愉快地用与DFT几乎完全相同的方法做IDFT了。
不过这个时候你可能并不满意——这个
V
−
1
V^{-1}
V−1 仿佛是直接从天上掉下来的,我虽然知道它是对的,但我还想知道:它是怎么推导出来的?
笑死,oi要什么推导过程,背板子就行了
不过这个看起来就挺让人发怵的玩意显然不在我当年学oi的时候的能力范围之内,于是直到前几天的数学分析习题课,我才知道了这东西的推导方法——几个爆算的推法和一个降维打击的推法,接下来一一展开介绍。
注:如果不想看繁琐的推导,可以直接跳到文章的最后一部分。
方法1:克莱姆法则
warning:流畅阅读这一部分内容需要一些线性代数前置知识。
warning2:这一部分推导过于繁琐。
这是习题课助教讲的第一个方法(大概讲了半个小时,可见有多么繁琐……)。
我们知道,线性方程组
A
x
=
b
Ax=b
Ax=b (其中
A
A
A 是
n
n
n 阶矩阵,
x
,
b
x,b
x,b是两个
n
n
n 维向量)有一个通用解法:
x
i
=
∣
A
i
∣
∣
A
∣
x_i=\frac{|A_i|}{|A|}
xi=∣A∣∣Ai∣ ,其中
A
i
A_i
Ai 是用向量
b
b
b 替换
A
A
A 中第
i
i
i 列得到的矩阵。这一方法称为克莱姆法则,证明在此略去。
我们直接利用这一方法求解,得
a
i
=
∣
V
i
∣
∣
V
∣
a_i=\frac{|V_i|}{|V|}
ai=∣V∣∣Vi∣。首先我们要求出
V
V
V 的行列式,容易看出
V
V
V 是范德蒙矩阵,因此有:
∣
V
∣
=
∏
0
≤
i
<
j
<
n
(
ω
j
−
ω
i
)
|V| = \prod_{0\leq i \lt j \lt n} (\omega^j - \omega^i)
∣V∣=∏0≤i<j<n(ωj−ωi)。这个式子不是那么容易化简,就先这么放着。
然后要求出
∣
V
i
∣
|V_i|
∣Vi∣ ,也就是:
d
e
t
(
1
1
⋯
1
b
0
1
⋯
1
1
ω
⋯
ω
i
−
1
b
1
ω
i
+
1
⋯
ω
n
−
1
1
ω
2
⋯
ω
2
(
i
−
1
)
b
2
ω
2
(
i
+
1
)
⋯
ω
2
(
n
−
1
)
⋮
⋮
⋮
⋮
⋮
⋮
1
ω
n
−
1
⋯
ω
(
n
−
1
)
(
i
−
1
)
b
n
−
1
ω
(
n
−
1
)
(
i
+
1
)
⋯
ω
(
n
−
1
)
2
)
det\left( \begin{matrix} 1 & 1 & \cdots & 1 & b_0 & 1 & \cdots & 1\\ 1 & \omega & \cdots & \omega^{i-1} & b_1 & \omega^{i+1} & \cdots & \omega^{n-1}\\ 1 & \omega^2 & \cdots & \omega^{2(i-1)} & b_2 & \omega^{2(i+1)} & \cdots & \omega^{2(n-1)}\\ \vdots & \vdots & & \vdots & \vdots & \vdots & & \vdots\\ 1 & \omega^{n-1} & \cdots & \omega^{(n-1)(i-1)} & b_{n-1} & \omega^{(n-1)(i+1)} & \cdots & \omega^{(n-1)^2}\\ \end{matrix} \right)
det⎝⎜⎜⎜⎜⎜⎛111⋮11ωω2⋮ωn−1⋯⋯⋯⋯1ωi−1ω2(i−1)⋮ω(n−1)(i−1)b0b1b2⋮bn−11ωi+1ω2(i+1)⋮ω(n−1)(i+1)⋯⋯⋯⋯1ωn−1ω2(n−1)⋮ω(n−1)2⎠⎟⎟⎟⎟⎟⎞
我们对被替换的列进行展开,得
∣
V
i
∣
=
∑
j
=
0
n
−
1
(
−
1
)
i
+
j
b
j
∣
V
j
i
∣
|V_i|=\sum_{j=0}^{n-1} (-1)^{i+j} b_j|V_{ji}|
∣Vi∣=∑j=0n−1(−1)i+jbj∣Vji∣ 。其中
V
j
i
V_{ji}
Vji 表示把矩阵
V
V
V 去掉第
j
j
j 行第
i
i
i 列得到的矩阵。
为了求出
∣
V
j
i
∣
|V_{ji}|
∣Vji∣,我们还需要求这样一类矩阵的行列式:
d
e
t
(
1
1
1
⋯
1
x
1
x
2
x
3
⋯
x
n
x
1
x
x
2
2
x
3
2
⋯
x
n
2
⋮
⋮
⋮
⋮
x
1
k
−
1
x
2
k
−
1
x
3
k
−
1
⋯
x
n
k
−
1
x
1
k
+
1
x
2
k
+
1
x
3
k
+
1
⋯
x
n
k
+
1
⋮
⋮
⋮
⋮
x
1
n
x
2
n
x
3
n
⋯
x
n
n
)
det\left( \begin{matrix} 1 & 1 & 1 & \cdots & 1\\ x_1 & x_2 & x_3 & \cdots & x_n\\ x_1^x & x_2^2 & x_3^2 & \cdots & x_n^2\\ \vdots & \vdots & \vdots & & \vdots\\ x_1^{k-1} & x_2^{k-1} & x_3^{k-1} & \cdots & x_n^{k-1}\\ x_1^{k+1} & x_2^{k+1} & x_3^{k+1} & \cdots & x_n^{k+1}\\ \vdots & \vdots & \vdots & & \vdots\\ x_1^n & x_2^n & x_3^n & \cdots & x_n^n\\ \end{matrix} \right)
det⎝⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎛1x1x1x⋮x1k−1x1k+1⋮x1n1x2x22⋮x2k−1x2k+1⋮x2n1x3x32⋮x3k−1x3k+1⋮x3n⋯⋯⋯⋯⋯⋯1xnxn2⋮xnk−1xnk+1⋮xnn⎠⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎞
这类矩阵长得很像范德蒙矩阵,我们在此不加证明(因为推导实在繁琐)地给出其行列式为
e
n
−
k
(
x
1
,
.
.
.
,
x
n
)
∏
(
1
≤
i
<
j
≤
n
)
(
x
j
−
x
i
)
e_{n-k}(x_1,...,x_n) \prod_{(1\leq i \lt j \leq n)} (x_j - x_i)
en−k(x1,...,xn)∏(1≤i<j≤n)(xj−xi) 。其中
e
k
(
x
1
,
.
.
.
.
,
x
n
)
e_k(x_1,....,x_n)
ek(x1,....,xn) 表示
x
1
,
.
.
.
,
x
n
x_1,...,x_n
x1,...,xn的第
k
k
k 个对称多项式,即:
e
k
(
x
1
,
.
.
.
.
,
x
n
)
=
∑
1
≤
i
1
<
i
2
<
⋯
<
i
k
≤
n
x
i
1
x
i
2
⋯
x
i
k
e_k(x_1,....,x_n)=\sum_{1\leq i_1\lt i_2 \lt \cdots \lt i_k \leq n} x_{i_1}x_{i_2}\cdots x_{i_k}
ek(x1,....,xn)=∑1≤i1<i2<⋯<ik≤nxi1xi2⋯xik。
然后带入化简,得:
a
i
=
∑
j
=
0
n
−
1
(
−
1
)
i
+
j
b
j
e
n
−
1
−
j
(
1
,
ω
,
⋯
,
ω
i
−
1
,
ω
i
+
1
,
⋯
,
ω
n
−
1
)
∏
0
≤
k
<
l
<
n
,
k
,
l
≠
i
(
ω
l
−
ω
k
)
∏
0
≤
k
<
l
<
n
(
ω
l
−
ω
k
)
a_i=\sum_{j=0}^{n-1} \frac{(-1)^{i+j}b_je_{n-1-j}(1,\omega,\cdots,\omega^{i-1},\omega^{i+1},\cdots,\omega^{n-1})\prod_{0\leq k \lt l\lt n,\ k,l\neq i}(\omega^l - \omega^k)}{\prod_{0\leq k \lt l\lt n}(\omega^l - \omega^k)}
ai=∑j=0n−1∏0≤k<l<n(ωl−ωk)(−1)i+jbjen−1−j(1,ω,⋯,ωi−1,ωi+1,⋯,ωn−1)∏0≤k<l<n, k,l=i(ωl−ωk)
方便起见,设
a
i
=
∑
j
=
0
n
−
1
c
i
j
b
j
a_i=\sum_{j=0}^{n-1} c_{ij}b_j
ai=∑j=0n−1cijbj,再记
e
i
,
j
=
e
i
(
1
,
ω
,
⋯
,
ω
j
−
1
,
ω
j
+
1
,
⋯
,
ω
n
−
1
)
e_{i,j}=e_i(1,\omega,\cdots,\omega^{j-1},\omega^{j+1},\cdots,\omega^{n-1})
ei,j=ei(1,ω,⋯,ωj−1,ωj+1,⋯,ωn−1) ,则有:
c
i
j
=
(
−
1
)
i
+
j
e
n
−
1
−
j
,
i
∏
i
<
k
<
n
(
ω
k
−
ω
i
)
∏
0
≤
k
<
i
(
ω
i
−
ω
k
)
c_{ij}=\frac{(-1)^{i+j}e_{n-1-j,i}}{\prod_{i\lt k \lt n} (\omega^k - \omega^i)\prod_{0 \leq k \lt i} (\omega^i - \omega^k)}
cij=∏i<k<n(ωk−ωi)∏0≤k<i(ωi−ωk)(−1)i+jen−1−j,i
=
(
−
1
)
n
+
j
+
1
e
n
−
1
−
j
,
i
∏
0
≤
k
<
n
,
k
≠
i
(
ω
i
−
ω
k
)
=\frac{(-1)^{n+j+1}e_{n-1-j,i}}{\prod_{0\leq k \lt n,\ k\neq i} (\omega^i - \omega^k)}
=∏0≤k<n, k=i(ωi−ωk)(−1)n+j+1en−1−j,i
接下来我们需要计算
e
i
,
j
e_{i,j}
ei,j。由
1
,
ω
,
⋯
,
w
n
−
1
1,\omega,\cdots,w^{n-1}
1,ω,⋯,wn−1 是方程
x
n
−
1
=
0
x^n - 1 = 0
xn−1=0 的根可得
x
n
−
1
=
(
x
−
1
)
(
x
−
ω
)
⋯
(
x
−
ω
n
−
1
)
x^n - 1 = (x - 1)(x - \omega)\cdots(x-\omega^{n-1})
xn−1=(x−1)(x−ω)⋯(x−ωn−1) ,对等式右侧直接展开并比较系数,有:
e
0
(
1
,
ω
,
.
.
.
,
ω
n
−
1
)
=
1
e_0(1,\omega,...,\omega^{n-1})=1
e0(1,ω,...,ωn−1)=1
e
n
(
1
,
ω
,
.
.
.
,
ω
n
−
1
)
=
(
−
1
)
n
+
1
e_n(1,\omega,...,\omega^{n-1})=(-1)^{n+1}
en(1,ω,...,ωn−1)=(−1)n+1
e
i
(
1
,
ω
,
.
.
.
,
ω
n
−
1
)
=
0
(
1
≤
i
<
n
)
e_i(1,\omega,...,\omega^{n-1})=0\ (1 \leq i \lt n)
ei(1,ω,...,ωn−1)=0 (1≤i<n)
再将
e
i
(
1
,
ω
,
.
.
.
,
ω
n
−
1
)
e_i(1,\omega,...,\omega^{n-1})
ei(1,ω,...,ωn−1) 中的项分为含有
ω
j
\omega^j
ωj 和不含有
ω
j
\omega^j
ωj 分开考虑,可得:
e
i
(
1
,
ω
,
.
.
.
,
ω
n
−
1
)
=
e
i
,
j
+
ω
j
e
i
−
1
,
j
e_i(1,\omega,...,\omega^{n-1})=e_{i,j}+\omega^j e_{i-1,j}
ei(1,ω,...,ωn−1)=ei,j+ωjei−1,j。
再结合
e
0
,
j
=
1
e_{0,j}=1
e0,j=1,可以求出
e
i
,
j
=
(
−
1
)
i
ω
i
j
e_{i,j}=(-1)^i\omega^{ij}
ei,j=(−1)iωij。
代回
c
i
j
c_{ij}
cij ,有:
c
i
j
=
(
−
1
)
n
+
j
+
1
(
−
1
)
n
−
j
−
1
ω
i
(
n
−
j
−
1
)
∏
0
≤
k
<
n
,
k
≠
i
(
ω
i
−
ω
k
)
=
ω
−
i
(
j
+
1
)
∏
0
≤
k
<
n
,
k
≠
i
(
ω
i
−
ω
k
)
c_{ij}=\frac{(-1)^{n+j+1}(-1)^{n-j-1}\omega^{i(n-j-1)}}{\prod_{0\leq k \lt n,\ k\neq i} (\omega^i - \omega^k)}=\frac{\omega^{-i(j+1)}}{\prod_{0\leq k \lt n,\ k\neq i} (\omega^i - \omega^k)}
cij=∏0≤k<n, k=i(ωi−ωk)(−1)n+j+1(−1)n−j−1ωi(n−j−1)=∏0≤k<n, k=i(ωi−ωk)ω−i(j+1)
最后只需要求出
∏
0
≤
k
<
n
,
k
≠
i
(
ω
i
−
ω
k
)
{\prod_{0\leq k \lt n,\ k\neq i} (\omega^i - \omega^k)}
∏0≤k<n, k=i(ωi−ωk) 即可。我们注意到该式其实就是多项式
∏
0
≤
k
<
n
,
k
≠
i
(
x
−
ω
k
)
{\prod_{0\leq k \lt n,\ k\neq i} (x - \omega^k)}
∏0≤k<n, k=i(x−ωk) 在
x
=
ω
i
x=\omega^i
x=ωi 处的取值,而由
x
n
−
1
=
∏
i
=
0
n
−
1
(
x
−
ω
i
)
x^n - 1 = \prod_{i=0}^{n-1}(x - \omega^i)
xn−1=∏i=0n−1(x−ωi) 我们可以立即得到
∏
0
≤
k
<
n
,
k
≠
i
(
x
−
ω
k
)
=
x
n
−
1
x
−
ω
i
{\prod_{0\leq k \lt n,\ k\neq i} (x - \omega^k)} = \frac{x^n-1}{x-\omega^i}
∏0≤k<n, k=i(x−ωk)=x−ωixn−1,直接做多项式除法就可以得到
x
n
−
1
x
−
ω
i
=
x
n
−
1
+
ω
i
x
n
−
2
+
⋯
+
ω
i
(
n
−
1
)
\frac{x^n-1}{x-\omega^i} = x^{n-1} + \omega^ix^{n-2}+\cdots+\omega^{i(n-1)}
x−ωixn−1=xn−1+ωixn−2+⋯+ωi(n−1),代入
x
=
ω
i
x=\omega^i
x=ωi,每一项的值都是
ω
i
(
n
−
1
)
\omega^{i(n-1)}
ωi(n−1) ,故结果为
n
ω
i
(
n
−
1
)
n\omega^{i(n-1)}
nωi(n−1)。
最后代回
c
i
j
c_{ij}
cij ,得
c
i
j
=
ω
−
i
(
j
+
1
)
n
ω
i
(
n
−
1
)
=
ω
−
i
j
n
c_{ij}=\frac{\omega^{-i(j+1)}}{n\omega^{i(n-1)}}=\frac{\omega^{-ij}}{n}
cij=nωi(n−1)ω−i(j+1)=nω−ij。
最后,设
n
n
n 阶矩阵
C
=
(
c
11
⋯
c
1
n
⋮
⋱
⋮
c
n
1
⋯
c
n
n
)
C=\left( \begin{matrix} c_{11} & \cdots & c_{1n}\\ \vdots & \ddots & \vdots\\ c_{n1} & \cdots & c_{nn} \end{matrix} \right)
C=⎝⎜⎛c11⋮cn1⋯⋱⋯c1n⋮cnn⎠⎟⎞ ,由
a
i
=
∑
j
=
0
n
−
1
c
i
j
b
j
a_i=\sum_{j=0}^{n-1} c_{ij}b_j
ai=∑j=0n−1cijbj ,可以立刻得到
α
=
C
β
\alpha=C\beta
α=Cβ,因此矩阵
C
C
C 就是我们要求的
V
−
1
V^{-1}
V−1 。
如果你存活到了现在,恭喜你数理基础很强。
这推导当然太繁琐了(而且里面很多定理的证明其实都跳过了),于是我们有接下来几种相对简单的方法:
方法2:范德蒙矩阵求逆
我记得习题课上助教提了一嘴这个算法,不过没有细算下去,因为你很快就会发现这和上面的爆算本质上没有区别。
对于一般的范德蒙矩阵
A
=
(
1
1
⋯
1
x
1
x
2
⋯
x
n
x
1
2
x
2
2
⋯
x
n
2
⋮
⋮
⋱
⋮
x
1
n
−
1
x
2
n
−
1
⋯
x
n
n
−
1
)
A=\left( \begin{matrix} 1 & 1& \cdots & 1\\ x_1 & x_2 & \cdots & x_n\\ x_1^2 & x_2^2 & \cdots & x_n^2\\ \vdots & \vdots & \ddots & \vdots\\ x_1^{n-1} & x_2^{n-1} & \cdots & x_n^{n-1} \end{matrix} \right)
A=⎝⎜⎜⎜⎜⎜⎛1x1x12⋮x1n−11x2x22⋮x2n−1⋯⋯⋯⋱⋯1xnxn2⋮xnn−1⎠⎟⎟⎟⎟⎟⎞,我们其实有通用的求逆公式:设
C
=
A
−
1
=
(
c
11
⋯
c
1
n
⋮
⋱
⋮
c
n
1
⋯
c
n
n
)
C=A^{-1}=\left( \begin{matrix} c_{11} & \cdots & c_{1n}\\ \vdots & \ddots & \vdots\\ c_{n1} & \cdots & c_{nn} \end{matrix} \right)
C=A−1=⎝⎜⎛c11⋮cn1⋯⋱⋯c1n⋮cnn⎠⎟⎞,则
c
i
j
=
(
−
1
)
j
−
1
e
n
−
j
(
x
1
,
x
2
,
⋯
,
x
i
−
1
,
x
i
+
1
,
⋯
,
x
n
)
∏
1
≤
k
≤
n
,
k
≠
i
(
x
k
−
x
i
)
c_{ij}=\frac{(-1)^{j-1}e_{n-j}(x_1,x_2,\cdots,x_{i-1},x_{i+1},\cdots,x_n)}{\prod_{1\leq k\leq n,k\neq i} (x_k-x_i)}
cij=∏1≤k≤n,k=i(xk−xi)(−1)j−1en−j(x1,x2,⋯,xi−1,xi+1,⋯,xn)。
然后你会发现接下来干的事跟方法1是一样的,无非是省下了一开始通过克莱姆法则得到
c
i
j
c_{ij}
cij 的过程(当然你也可以看作方法1是用克莱姆法则强行推了一遍范德蒙矩阵求逆公式)。
方法3:拉格朗日插值
这是我在习题课上yy出来的推法。
思路很简单,我们既然已经有了
f
(
1
)
,
f
(
ω
)
,
f
(
ω
2
)
,
⋯
,
f
(
ω
n
−
1
)
f(1),f(\omega),f(\omega^2),\cdots,f(\omega^{n-1})
f(1),f(ω),f(ω2),⋯,f(ωn−1),本质上就是有了
f
f
f 这个
n
−
1
n-1
n−1 次多项式在
n
n
n 个点处的点值。我们直接套用拉格朗日插值公式,就可以得到:
f
(
x
)
=
∑
i
=
0
n
−
1
b
i
∏
j
≠
i
(
x
−
ω
j
)
∏
j
≠
i
(
ω
i
−
ω
j
)
f(x)=\sum_{i=0}^{n-1}b_i\frac{\prod_{j\neq i}(x-\omega^j)}{\prod_{j\neq i}(\omega^i - \omega^j)}
f(x)=∑i=0n−1bi∏j=i(ωi−ωj)∏j=i(x−ωj) 。
其中这个分子和分母我们都已经算出来了:
∏
j
≠
i
(
x
−
ω
j
)
=
x
n
−
1
+
ω
i
x
n
−
2
+
⋯
+
ω
i
(
n
−
1
)
=
∑
j
=
0
n
−
1
ω
i
(
n
−
j
−
1
)
x
j
,
∏
j
≠
i
(
ω
i
−
ω
j
)
=
n
ω
i
(
n
−
1
)
\prod_{j\neq i}(x-\omega^j)=x^{n-1}+\omega^ix^{n-2}+\cdots+\omega^{i(n-1)}=\sum_{j=0}^{n-1}\omega^{i(n-j-1)}x^j,\ \prod_{j\neq i}(\omega^i - \omega^j)=n\omega^{i(n-1)}
∏j=i(x−ωj)=xn−1+ωixn−2+⋯+ωi(n−1)=∑j=0n−1ωi(n−j−1)xj, ∏j=i(ωi−ωj)=nωi(n−1) 。
于是可以立刻得到
f
(
x
)
=
∑
i
=
0
n
−
1
b
i
∑
j
=
0
n
−
1
ω
−
i
j
n
x
j
=
∑
i
=
0
n
−
1
(
1
n
∑
j
=
0
n
−
1
ω
−
i
j
b
j
)
x
i
f(x)=\sum_{i=0}^{n-1}b_i\sum_{j=0}^{n-1}\frac{\omega^{-ij}}{n}x^j=\sum_{i=0}^{n-1}(\frac{1}{n}\sum_{j=0}^{n-1}\omega^{-ij}b_j)x^i
f(x)=∑i=0n−1bi∑j=0n−1nω−ijxj=∑i=0n−1(n1∑j=0n−1ω−ijbj)xi ,换句话说
a
i
=
1
n
∑
j
=
0
n
−
1
ω
−
i
j
b
j
a_i=\frac{1}{n}\sum_{j=0}^{n-1}\omega^{-ij}b_j
ai=n1∑j=0n−1ω−ijbj 。
这种方法就简单多了,因为不需要繁琐的矩阵运算,也不需要计算对称多项式的值,只需要处理一下
x
n
−
1
x
−
ω
i
\frac{x^n-1}{x-\omega^i}
x−ωixn−1 就可以了。
方法4:复内积
这是习题课上助教最后讲的“降维打击”,我们可以看到原本极其繁琐的计算是怎样通过引入内积之后瞬间解决的。
设
X
=
(
x
1
,
x
2
,
⋯
,
x
n
)
,
Y
=
(
y
1
,
y
2
,
⋯
,
y
n
)
X=(x_1,x_2,\cdots,x_n),\ Y=(y_1,y_2,\cdots,y_n)
X=(x1,x2,⋯,xn), Y=(y1,y2,⋯,yn) 是两个
n
n
n 维复向量,定义它们的复内积为
<
X
,
Y
>
=
∑
i
=
1
n
x
i
y
i
‾
<X,\ Y>\ =\sum_{i=1}^n x_i\overline{y_i}
<X, Y> =∑i=1nxiyi。
将最初的矩阵
V
V
V看成
n
n
n 个列向量并排的形式:
V
=
(
γ
0
γ
1
⋯
γ
n
−
1
)
V=(\gamma_0\ \gamma_1\ \cdots \ \gamma_{n-1})
V=(γ0 γ1 ⋯ γn−1) ,其中
γ
i
=
(
1
ω
i
ω
2
i
⋮
ω
i
(
n
−
1
)
)
\gamma_i=\left( \begin{matrix} 1\\ \omega^i\\ \omega^{2i}\\ \vdots\\ \omega^{i(n-1)} \end{matrix} \right)
γi=⎝⎜⎜⎜⎜⎜⎛1ωiω2i⋮ωi(n−1)⎠⎟⎟⎟⎟⎟⎞。
则有
β
=
∑
i
=
0
n
−
1
a
i
γ
i
\beta = \sum_{i=0}^{n-1}a_i\gamma_i
β=∑i=0n−1aiγi 。
接下来一步是神来之笔:对于每个
i
=
0
,
⋯
,
n
−
1
i=0,\cdots,n-1
i=0,⋯,n−1,等号两边同时与
γ
i
\gamma_i
γi做复内积:
<
β
,
γ
i
>
=
<
∑
j
=
0
n
−
1
a
j
γ
j
,
γ
i
>
=
∑
j
=
0
n
−
1
a
i
<
γ
j
,
γ
i
>
<\beta,\gamma_i>\ =\ <\sum_{j=0}^{n-1}a_j\gamma_j,\gamma_i>\ =\sum_{j=0}^{n-1}a_i<\gamma_j,\gamma_i>
<β,γi> = <∑j=0n−1ajγj,γi> =∑j=0n−1ai<γj,γi>。
其中,
<
γ
i
,
γ
j
>
=
∑
k
=
0
n
−
1
ω
i
k
ω
j
k
‾
=
∑
k
=
0
n
−
1
ω
(
i
−
j
)
k
<\gamma_i,\gamma_j>\ =\sum_{k=0}^{n-1}\omega^{ik}\overline{\omega^{jk}}=\sum_{k=0}^{n-1}\omega^{(i-j)k}
<γi,γj> =∑k=0n−1ωikωjk=∑k=0n−1ω(i−j)k。
当
i
=
j
i=j
i=j 时,显然有
<
γ
i
,
γ
j
>
=
n
<\gamma_i,\gamma_j>\ =n
<γi,γj> =n,否则由对称性容易看出
<
γ
i
,
γ
j
>
=
0
<\gamma_i,\gamma_j>\ =0
<γi,γj> =0。
因此,上述
<
β
,
γ
i
>
<\beta,\gamma_i>
<β,γi> 的等式右边就只剩
<
γ
i
,
γ
i
>
<\gamma_i,\gamma_i>
<γi,γi> 一项了,于是我们得到
<
β
,
γ
i
>
=
n
a
i
<\beta,\gamma_i>\ = na_i
<β,γi> =nai。整理得
a
i
=
1
n
∑
j
=
0
n
−
1
ω
−
i
j
b
j
a_i=\frac{1}{n}\sum_{j=0}^{n-1}\omega^{-ij}b_j
ai=n1∑j=0n−1ω−ijbj 。
这当然是极好的!我们几乎没有用任何线性代数知识,没有任何文中没有证明的前置知识,没有任何繁琐的计算和推导,仅引入了一个数学工具就成功解决了问题!
什么,你说你还是看不懂?那我建议你先把学FFT的事一放,先去巩固一下你的数理基础吧