Linear Algebra-酉矩阵与快速FFT-27

一、知识概要

在上一节中，我们提到过复数矩阵，但并没有深入讨论它的运算与性质。本节我们就来补上这一课，一起学习复数矩阵的基本运算和特征，然后重点认识一个非常重要的复数矩阵——傅里叶矩阵。最后，我们会介绍大名鼎鼎的快速傅里叶变换（FFT），看看它是如何显著减少运算量的。

二、复数矩阵

我们从最基础的复向量开始说起。如果向量中的分量出现了虚数，我们就不能再用实向量的公式来计算长度和内积了。举个例子，比如向量，如果按照实向量的方式计算：，算出来长度竟然是0，但实际上我们知道这个向量的模长应该是，这就产生了矛盾。所以，我们需要重新给出复向量的长度与内积计算方法。

对于复数向量，向量长度的计算公式为，同样，两个复数向量的内积定义为。这里的上标表示先取共轭再转置，也叫作共轭转置。

接下来我们聊聊对称阵。在实矩阵中，若，我们就称为对称阵。但在复矩阵中，这个定义就不再适用了。结合上面复向量内积的处理方法，我们会发现，在复矩阵中取共轭与取转置常常是同时进行的。因此，复对称矩阵的定义是：若，则称为对称阵，也叫埃尔米特矩阵，记作。

比如矩阵，主对角线上都是实数，沿主对角线对称的两个元素正好共轭，这就是一个典型的埃尔米特矩阵。埃尔米特矩阵有两个非常好的性质：特征值都是实数，且不同特征值对应的特征向量正交。

同理，我们把特征向量正交的概念也推广到复矩阵中。若一组向量满足，我们就称它们是一组标准正交基。把这组标准正交基按列排列构成矩阵，即，那么我们就有：

满足这个性质的矩阵就称为酉矩阵。这相当于实矩阵中的正交矩阵，只不过推广到了复数域。

三、傅里叶矩阵与快速傅里叶变换

3.1 傅里叶矩阵

傅里叶矩阵是一个非常典型的酉矩阵，它的形式是这样的：

其中，且。

计算的时候，我们一般都用这种指数形式，而不是写成的形式。从几何意义上看，表示单位圆上的一个点，就是把单位圆等分的份数， ， ， ， ，正好对应这些等分点。

我们举个最简单的例子，假设，那么，对应的四阶傅里叶矩阵就是：

可以验证一下，这个矩阵的各个列向量确实是正交的，满足。如果我们把乘以得到新的，那么新矩阵就满足，这就是一个标准的酉矩阵了，它的逆矩阵就是它的共轭转置。

3.2 快速傅里叶变换

快速傅里叶变换（FFT）的核心思想非常巧妙：大尺寸的傅里叶矩阵和更小尺寸的傅里叶矩阵之间存在着分解关系，利用这种分解关系可以大幅减少计算量。我们以和为例来说明：

直接代入公式，我们会发现，中的正好等于，这就建立了两者之间的联系。我们可以构造出如下的矩阵分解公式：

我们来解释一下这个式子各个部分的作用：

  • 置换矩阵：作用是将矩阵的奇偶列（或奇偶行）分开排列，这是为了分解计算做准备；
  • 对角矩阵：的对角元素是，其中，用来修正不同位置的系数。

举个小例子，在4维情况下，置换矩阵，它把原来的第2行换到了第3行的位置，实现了奇偶分离。而，整个的作用就是把两个小块的组合成完整的。

现在我们来算算计算量的变化。如果直接用乘以一个列向量，因为是64×64的矩阵，所以原本需要次运算。而经过上面的分解，计算量变成了多少呢？两个的运算量是次，再加上修正项中的32阶对角矩阵，需要额外32次运算。这里我们注意到，虽然式子中出现了和，但它们只是符号相反，具体数值只需要计算一次，改变符号就能得到另一个结果，所以只算一次就够了。这样总计算量就是次，已经比原来少了很多。

我们还可以继续分解下去，把每个再分解成两个加上修正项，计算量会进一步减小为次。就这样一直分解下去，一共要做次分解，最后傅里叶矩阵就变成了单位矩阵，只需要计算修正项了。

最终，对于，总的计算量是次。这个结果太惊艳了！推广到一般的n阶傅里叶矩阵，快速傅里叶变换的总运算量为，和原来的次运算相比，减少得非常显著。我们画张图感受一下这个变化：

当n很大的时候，比如n=1024，原来需要100多万次运算，FFT只需要大约5000次，差距一目了然。

3.3 验证用的 Matlab 代码

我们可以用Matlab来验证一下上面的分解公式是否正确，下面是完整的代码：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
% 计算傅里叶矩阵
% 计算F64
n = 64;
w = exp(1j*2*pi/n);
F64 = zeros(n);
for i = 1:n
    for j = 1:n
        F64(i,j) = w^((i - 1)*(j - 1));
    end
end

% 计算F32
m = 32;
w32 = exp(1j*2*pi/m);
F32 = zeros(m);
for i = 1:m
    for j = 1:m
        F32(i,j) = w32^((i - 1)*(j - 1));
    end
end

% 构造置换矩阵P
P = eye(n);
P = P([1:2:n, 2:2:n], :);

% 构造对角矩阵D
D = diag(w.^(0:m - 1));

% 计算等式右边的矩阵
right_hand_side = [eye(m), D; eye(m), -D]*[F32, zeros(m); zeros(m), F32]*P;

% 验证等式是否成立
error_norm = norm(F64 - right_hand_side);
if error_norm < 1e-10
    disp('理论验证通过，两个矩阵在数值上相等。');
else
    disp('理论验证未通过，两个矩阵在数值上不相等。');
end

运行这段代码，我们就能验证分解公式的正确性了。

三、学习感悟

本节主要研究了复矩阵、复向量的计算以及性质，拓展了向量内积和向量长度的计算方法，还探讨了一种特殊的复矩阵——傅里叶矩阵。最后讲解了本章重点内容快速傅里叶变换，这部分是难点，尤其是对这一转换的理解有一定难度。个人理解是转换矩阵将向量的奇偶行分离，从而减小计算量，而的作用是实现两个傅里叶矩阵之间的转换。快速傅里叶变换（FFT）非常重要，其应用较为复杂，想要彻底理解还需要进一步深入学习这部分知识。

四、学习总结

  1. 埃尔米特矩阵是复数域的对称矩阵，；
  2. 埃尔米特矩阵的特征值都是实数，且特征向量正交；
  3. 酉矩阵是复数域的标准正交矩阵，；