傅里叶变换的差分矩阵

                     

贡献者： 零穹; addis

　　 差分矩阵 $D$ 和数据点矢量 $u$ 矩阵乘法 $Du$ 给出数据点处的近似导数。本文将给出基于傅里叶变换给出的差分矩阵。

1. 傅里叶变换回顾

连续区间上的傅里叶变换

　　 由正交函数系中的例 3 可知，函数 $u(x),x\in \mathbb{R}$ 的傅里叶变换定义为

而从 $\hat{u}$ 的逆傅里叶变换可以重构 $u$： 这被称为傅里叶合成（Fourier synthesis），变量 $x$ 称为物理变量（physical variable），$k$ 称为傅里叶变量（Fourier variable）或波数（wavenumber）。

离散点上的傅里叶变换

　　 无界情形

　　 当限定 $x\in h\mathbb{Z}$ 时，即此时 $x$ 只能取离散点，那么只要 $k_1-k_2=\frac{2\pi }{h}$，就有 ${\mathrm{e}}^{\mathrm{i}{k}_{1}x}={\mathrm{e}}^{\mathrm{i}{k}_{2}x}$ 对所有的 $x$ 成立。这使得对任意复指数 ${\mathrm{e}}^{\mathrm{i}kx}$，有无穷多个复指数和它在 $h\mathbb{Z}$ 上具有相同的值。这一现象被称为混叠。为了保持式 1 中傅里叶变换 $\hat{u}(k)$ 的单值性，需要将 $k$ 限制在长为 $\frac{2\pi }{h}$ 的区间上。因此，我们可以选择 $k$ 的取值区间为 $[-\pi /h,\pi /h]$。注意到 $x$ 是物理变量，$k$ 是傅里叶变量，因此我们可以将这些结果总结为下面的图示：

因此，在无界离散的情形，傅里叶变换及其逆变换式 1 和式 2 应当写为下面的级数和： 这被称为半离散傅里叶变换。

　　 周期性情形

　　 当限制 $u(x_i)$ 对 $x_i$ 具有周期 $2L$ 时，即 $u(x_i+2L)=u(x_i)$，那么为了 $u$ 的单值性，限定 $x_i\in [0,2L]$。设 $x_i$ 有 $N$ 个（本词条将始终假定 $N$ 为偶数），于是 $x_j\in \{h,\cdots ,Nh=2L\}$，其中 $h:=\frac{2L}{N}$。仍记 $x_j:=jh,j\in \{1,\cdots ,N\}$。和无界情形一样的理由，$k$ 应当限制在 $[-\pi /h,\pi /h]=[-\frac{N\pi }{2L},\frac{N\pi }{2L}]$ 上。此外由于周期性，应当认为点 $x+2L$ 和点 $x$ 是同一个。因此 ${\mathrm{e}}^{\mathrm{i}k(x+2L)}={\mathrm{e}}^{\mathrm{i}kx}$，从而

结合 $k\in [-\frac{N\pi }{2L},\frac{N\pi }{2L}]$，有 注意这里用到了 $-\frac{N}{2}\frac{\pi }{L},\frac{N}{2}\frac{\pi }{L}$ 对应相同的 ${\mathrm{e}}^{-\mathrm{i}kx}$。这些结果可以总结为下面的图示： 相应的，傅里叶变换和逆变换的公式写为下面的级数和， 这被称为离散傅里叶变换。

2. 插值函数

无界情形

　　 假若我们有定义在 $x_i\in h\mathbb{Z}$ 上的数据 $u_i=u(x_i)$，那么式 4 将会给我们插值函数：

　　 1.通过式 4 的第一式我们可以获得 $\hat{u}(k)$，它关于 $k$ 是可微的;

　　 2.而式 4 的第二式可以确定一个函数

显然 $p(x)$ 满足 $p(x_i)=u_i$，并且 $p(x)$ 是可微的。因此可以用 $p(x)$ 的导数来近似 $u(x)$ 的导数。$p(x)$ 给了我们一个插值公式。

　　 由式 9 ，只需令 $\hat{u}(k)$ 仅在 $[-\pi /h,\pi /h]$ 上不为零。那么式 9 的积分下限和上限就可以变到 $-\mathrm{\infty },\mathrm{\infty }$。因此由式 1 和式 2 可知 $p(x)$ 的傅里叶变换为

由于这一点，式 9 的 $p(x)$ 被称为是 $v$ 的有限带宽插值函数。这一术语暗示着 $p$ 的傅里叶变换仅在 $[-\pi /h,\pi /h]$ 上不为 0（仅在某一区域 $U$ 上不为 0 的函数往往称为具有支撑，$U$ 称为它的支撑集。若 $U$ 还是拓扑意义上的紧集，又被称为具有紧支撑）。

周期情形

　　 假设数据定义在 $x_i=ih,i\in \{1,\cdots N\},h=\frac{2L}{N}$ 上，$u_i=u(x_i)$。那么：

　　 1.通过式 8 的第一式我们可以获得 $\hat{u}(k)$.

　　 2.而式 8 的第二式可以确定一个函数（为了更对称的处理，并考虑到 ${\mathrm{e}}^{\mathrm{i}k\frac{\pi }{L}x}$ 在 $k=\frac{N}{2},-\frac{N}{2}$ 相等）

其中 ${}^{\prime }$ 表示求和在 $-N/2,N/2$ 处需要乘 $1/2$。显然 $p(x)$ 满足 $p(x_i)=u_i$，并且 $p(x)$ 是可微的。因此可以用 $p(x)$ 的导数来近似 $u(x)$ 的导数。$p(x)$ 便是需要的插值公式。

3. 导数近似算法

无界情形

　　 有了式 9 的插值函数 $p(x)$，我们就可以进行导数近似了：

　　 1.给定定义在 $h\mathbb{Z}$ 上的函数 $u$，利用式 9 确定它的有限带宽函数 $p(x)$；

　　 2.令 $w_j=p^{\prime }(x_j)$，则 $w_j$ 就是 $u$ 在 $x_j$ 处的近似导数。

　　 当然，我们还有另一等价的方式来描述，它是在傅里叶空间中表述的：式 2 对 $x$ 微分得

由傅里叶正反变换的关系式 1 和式 2 可得 因此我们得到另一等价的导数近似程序：

　　 1.给定 $u$，利用式 4 的第一式计算它的离散傅里叶变换 $\hat{u}$；

　　 2.定义 $\hat{w}(k)=\mathrm{i}k\hat{u}(k)$;

　　 3.利用式 4 的第二式从 $\hat{w}$ 计算 $w$。这即得到在给定点的导数近似。

周期情形

　　 周期情形的程序为：

　　 1.给定定义在 $\{h,\cdots ,Nh=2L\}$（$x_i=ih,i=1,\cdots ,N$）上的函数 $u$，利用式 11 确定它的插值函数 $p(x)$；

　　 2.令 $w_j=p^{\prime }(x_j)$，则 $w_j$ 就是 $u$ 在 $x_j$ 处的近似导数。

4. 差分矩阵

无界情形

　　 由式 4 可知，若 $u(x)=\sum _j{a}_j{f}_j(x)$，$a_i\in \mathbb{C}$，那么

即傅里叶变换是线性的。而任意的离散点 $u(x_j)$ 都可以用 Kronecker Delta 函数（$x_j=jh$） 表示为 $u(x_j)=\sum _{m=-\mathrm{\infty }}^{\mathrm{\infty }}u(x_m)\delta (x_j-x_m)$。因此，若知道 $\delta (x-x_m)$ 的傅里叶变换 ${\hat{\delta }}_m(k)$，则式 4 的第一式就为 从而式 9 的 $p(x)$ 成为

　　 而

所以 其中 $\mathrm{sinc}(x):=\frac{\sin x}{x}$。从而 因为 $x\ne 0$ 时， 而 ${\mathrm{sinc}}^{\prime }0=0$，所以 $x\ne x_m$ 时 所以式 20 成为 该式表明，若已知数据点 $u_i=u(x_i),i\in \mathbb{Z}$，则其对应的导数可由上式求得。由子节 2 中差分矩阵的概念，即差分矩阵和数据点矢量的矩阵乘法给出数据点的近似导数，因此我们得到定义在 $h\mathbb{Z}$ 上对应的差分矩阵，其是无穷阶的（式 23 表明对应指标(j,m) 的矩阵元为 $\frac{(-1{)}^{j-m}}{(j-m)h}$，容易验证其是Toeplitz 矩阵）：

周期情形

　　 同样的，式 11 可以写为

而 所以

　　 由

　　 从而

因此，差分矩阵为