贡献者： 零穹
差分矩阵 $D$ 和数据点矢量 $u$ 矩阵乘法 $Du$ 给出数据点处的近似导数。本词条将给出基于傅里叶变换给出的差分矩阵。
1. 傅里叶变换回顾
连续区间上的傅里叶变换
由正交函数系中的例 3 可知，函数 $u(x),x\in\mathbb R$ 的傅里叶变换定义为
\begin{equation}
\hat u(k):=\int_{-\infty}^{\infty} e^{- \mathrm{i} kx}u(x) \,\mathrm{d}{x} ,\quad k\in\mathbb R.~
\end{equation}
而从 $\hat u$ 的逆傅里叶变换可以重构 $u$：
\begin{equation}
u(x)=\frac{1}{2\pi}\int_{-\infty}^{\infty} e^{ \mathrm{i} kx}\hat u(k) \,\mathrm{d}{k} .~
\end{equation}
这被称为
傅里叶合成（Fourier synthesis），变量 $x$ 称为
物理变量（physical variable），$k$ 称为
傅里叶变量（Fourier variable）或
波数（wavenumber）。
离散点上的傅里叶变换
无界情形
当限定 $x\in h\mathbb Z$ 时，即此时 $x$ 只能取离散点，那么只要 $k_1-k_2=\frac{2\pi}{h}$，就有 $e^{ \mathrm{i} k_1x}=e^{ \mathrm{i} k_2x}$ 对所有的 $x$ 成立。这使得对任意复指数 $e^{ \mathrm{i} kx}$，有无穷多个复指数和它在 $h\mathbb Z$ 上具有相同的值。这一现象被称为混叠。为了保持式 1 中傅里叶变换 $\hat u(k)$ 的单值性，需要将 $k$ 限制在长为 $\frac{2\pi}{h}$ 的区间上。因此，我们可以选择 $k$ 的取值区间为 $[-\pi/h,\pi/h]$。注意到 $x$ 是物理变量，$k$ 是傅里叶变量，因此我们可以将这些结果总结为下面的图示：
\begin{equation}
\begin{aligned}
&\text{物理空间:}& \text{离散},&\qquad\text{无界:} &x\in h\mathbb Z\\
&&\updownarrow&\qquad\updownarrow&\\
&\text{傅里叶空间:}& \text{有界},&\qquad\text{连续:} &k\in [-\pi/h,\pi/h]
\end{aligned}~
\end{equation}
因此，在无界离散的情形，傅里叶变换及其逆变换
式 1 和
式 2 应当写为下面的级数和：
\begin{equation}
\begin{aligned}
\hat u(k)&=h\sum_{i=-\infty}^{\infty}e^{- \mathrm{i} kx_i}u(x_i),\quad k\in[-\pi/h,\pi/h],\\
u(x_j)&=\frac{1}{2\pi}\int_{-\pi/h}^{\pi/h} e^{ \mathrm{i} kx_j}\hat u(k) \,\mathrm{d}{k} ,\quad x_j=jh,j\in \mathbb Z.
\end{aligned}~
\end{equation}
这被称为
半离散傅里叶变换。
周期性情形
当限制 $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$ 是同一个。因此 $e^{ \mathrm{i} k(x+2L)}=e^{ \mathrm{i} k x}$，从而
\begin{equation}
2Lk=2n\pi\Rightarrow k=n\frac{\pi}{L},\quad n\in\mathbb Z.~
\end{equation}
结合 $k\in[-\frac{N\pi}{2L},\frac{N\pi}{2L}]$，有
\begin{equation}
k\in \left\{ \left(-\frac{N}{2}+1 \right) \frac{\pi}{L}, \left(-\frac{N}{2}+2 \right) \frac{\pi}{L},\cdots,\frac{N}{2}\frac{\pi}{L} \right\} .~
\end{equation}
注意这里用到了 $-\frac{N}{2}\frac{\pi}{L},\frac{N}{2}\frac{\pi}{L}$ 对应相同的 $e^{- \mathrm{i} kx}$。这些结果可以总结为下面的图示：
\begin{equation}
\begin{aligned}
&\text{物理空间:}& \text{离散},&\qquad\text{有界:} \qquad x\in \{h,\cdots,2L\}\\
&&\updownarrow&\qquad\updownarrow&\\
&\text{傅里叶空间:}& \text{有界},&\qquad\text{离散:} \qquad k\in \left\{ \left(-\frac{N}{2}+1 \right) \frac{\pi}{L}, \left(-\frac{N}{2}+2 \right) \frac{\pi}{L},\cdots,\frac{N}{2}\frac{\pi}{L} \right\} .
\end{aligned}~
\end{equation}
相应的，傅里叶变换和逆变换的公式写为下面的级数和，
\begin{equation}
\begin{aligned}
\hat u(k)&=h\sum_{i=1}^{N}e^{- \mathrm{i} k\frac{\pi}{L}x_i}u(x_i),\quad k\in \left\{ \left(-\frac{N}{2}+1 \right) , \left(-\frac{N}{2}+2 \right) ,\cdots,\frac{N}{2} \right\} ,\\
u(x_j)&=\frac{1}{2L}\sum_{k=-\frac{N}{2}+1}^{\frac{N}{2}}e^{ \mathrm{i} k\frac{\pi}{L}x_j}\hat u(k),\quad j\in\{1,\cdots,N\}.
\end{aligned}~
\end{equation}
这被称为
离散傅里叶变换。
2. 插值函数
无界情形
假若我们有定义在 $x_i\in h\mathbb Z$ 上的数据 $u_i=u(x_i)$，那么式 4 将会给我们插值函数：
1.通过式 4 的第一式我们可以获得 $\hat u(k)$，它关于 $k$ 是可微的;
2.而式 4 的第二式可以确定一个函数
\begin{equation}
p(x)=\frac{1}{2\pi}\int_{-\pi/h}^{\pi/h} e^{ \mathrm{i} kx}\hat u(k) \,\mathrm{d}{k} ,\quad x\in \mathbb R.~
\end{equation}
显然 $p(x)$ 满足 $p(x_i)=u_i$，并且 $p(x)$ 是可微的。因此可以用 $p(x)$ 的导数来近似 $u(x)$ 的导数。$p(x)$ 给了我们一个插值公式。
由式 9 ，只需令 $\hat u(k)$ 仅在 $[-\pi/h,\pi/h]$ 上不为零。那么式 9 的积分下限和上限就可以变到 $-\infty,\infty$。因此由式 1 和式 2 可知 $p(x)$ 的傅里叶变换为
\begin{equation}
\hat p(k)=\left\{\begin{aligned}
&\hat u(k),&k\in[-\pi/h,\pi/h],\\
&0,&k\not\in[-\pi/h,\pi/h].
\end{aligned}\right.~
\end{equation}
由于这一点，
式 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 的第二式可以确定一个函数（为了更对称的处理，并考虑到 $e^{ \mathrm{i} k\frac{\pi}{L}x}$ 在 $k=\frac{N}{2},-\frac{N}{2}$ 相等）
\begin{equation}
p(x)=\frac{1}{2L}{\sum'}\limits_{k=-\frac{N}{2}}^{\frac{N}{2}}e^{ \mathrm{i} k\frac{\pi}{L}x}\hat u(k),~
\end{equation}
其中 $'$ 表示求和在 $-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'(x_j)$，则 $w_j$ 就是 $u$ 在 $x_j$ 处的近似导数。
当然，我们还有另一等价的方式来描述，它是在傅里叶空间中表述的：式 2 对 $x$ 微分得
\begin{equation}
u'(x)=\frac{1}{2\pi}\int_{-\infty}^{\infty} \mathrm{i} ke^{ \mathrm{i} kx}\hat u(k) \,\mathrm{d}{k} ,\quad x\in \mathbb R.~
\end{equation}
由傅里叶正反变换的关系
式 1 和
式 2 可得
\begin{equation}
\hat{u'}(k)=ik\hat u(k).~
\end{equation}
因此我们得到另一等价的导数近似程序：
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'(x_j)$，则 $w_j$ 就是 $u$ 在 $x_j$ 处的近似导数。
4. 差分矩阵
无界情形
由式 4 可知，若 $u(x)=\sum\limits_{j}a_jf_j(x)$，$a_i\in\mathbb C$，那么
\begin{equation}
\begin{aligned}
\hat u(k)&=h\sum_{i=-\infty}^{\infty}e^{- \mathrm{i} kx_i} \left(\sum\limits_{j}a_jf_j(x_i) \right) \\
&=\sum\limits_{j}a_j \left(h\sum_{i=-\infty}^{\infty}e^{- \mathrm{i} kx_i}f_j(x_i) \right) \\
&=\sum_{j}a_j\hat f_j(k).
\end{aligned}~
\end{equation}
即傅里叶变换是线性的。而任意的离散点 $u(x_j)$ 都可以用 Kronecker Delta 函数（$x_j=jh$）
\begin{equation}
\delta(x)=\left\{\begin{aligned}
&1,&x=0,\\
&0,&x\neq0.
\end{aligned}\right.~
\end{equation}
表示为 $u(x_j)=\sum\limits_{m=-\infty}^{\infty}u(x_m)\delta(x_j-x_m)$。因此，若知道 $\delta(x-x_m)$ 的傅里叶变换 $\hat\delta_m(k)$，则
式 4 的第一式就为
\begin{equation}
\hat u(k)=\sum\limits_{m=-\infty}^\infty u(x_m)\hat\delta_m(k).~
\end{equation}
从而
式 9 的 $p(x)$ 成为
\begin{equation}
\begin{aligned}
p(x)&=\frac{1}{2\pi}\int_{-\pi/h}^{\pi/h} e^{ \mathrm{i} kx}\hat u(k) \,\mathrm{d}{k} \\
&=\sum\limits_{m=-\infty}^\infty u(x_m)\frac{1}{2\pi}\int_{-\pi/h}^{\pi/h} e^{ \mathrm{i} kx}\hat\delta_m(k) \,\mathrm{d}{k} .
\end{aligned}~
\end{equation}
而
\begin{equation}
\begin{aligned}
\hat\delta_m(k)=h\sum_{i=-\infty}^{\infty}e^{- \mathrm{i} kx_i}\delta(x_i-x_m)=he^{- \mathrm{i} kx_m}.
\end{aligned}~
\end{equation}
所以
\begin{equation}
\begin{aligned}
p(x)&=\sum\limits_{m=-\infty}^\infty u(x_m)\frac{h}{2\pi}\int_{-\pi/h}^{\pi/h} e^{ \mathrm{i} k(x-x_m)} \,\mathrm{d}{k} \\
&=\sum\limits_{m=-\infty}^\infty u(x_m)\frac{h}{2\pi}\frac{e^{ \mathrm{i} \frac{\pi}{h}(x-x_m)}-e^{- \mathrm{i} \frac{\pi}{h}(x-x_m)}}{ \mathrm{i} (x-x_m)}\\
&=\sum\limits_{m=-\infty}^\infty u(x_m)\frac{h}{\pi}\frac{ \sin\left(\frac{\pi}{h}(x-x_m)\right) }{x-x_m}\\
&=\sum\limits_{m=-\infty}^\infty u(x_m) \operatorname{sinc} \left(\frac{\pi}{h}(x-x_m) \right) ,
\end{aligned}~
\end{equation}
其中 $ \operatorname{sinc} (x):=\frac{\sin x}{x}$。从而
\begin{equation}
\begin{aligned}
w_j&=p'(x_j)=\sum\limits_{m=-\infty}^\infty u(x_m) \operatorname{sinc} ' \left(\frac{\pi}{h}(x-x_m) \right) |_{x=x_j=jh}.
\end{aligned}~
\end{equation}
因为 $x\neq0$ 时，
\begin{equation}
\operatorname{sinc} '(x)= \left(\frac{\sin x}{x} \right) '=\frac{x\cos x-\sin x}{x^2},~
\end{equation}
而 $ \operatorname{sinc} '0=0$，所以 $x\neq x_m$ 时
\begin{equation}
\begin{aligned}
& \operatorname{sinc} ' \left(\frac{\pi}{h}(x-x_m) \right) |_{x=x_j=jh}=\frac{\pi}{h} \operatorname{sinc} '[(j-m)\pi]\\
&=\frac{(-1)^{j-m}}{(j-m)h}.
\end{aligned}~
\end{equation}
所以
式 20 成为
\begin{equation}
w_j=p'(x_j)=\sum\limits_{m=-\infty,m\neq j}^\infty \frac{(-1)^{j-m}}{(j-m)h}u(x_m).~
\end{equation}
该式表明，若已知数据点 $u_i=u(x_i),i\in\mathbb Z$，则其对应的导数可由上式求得。由
子节 2
中差分矩阵的概念，即差分矩阵和数据点矢量的矩阵乘法给出数据点的近似导数，因此我们得到定义在 $h\mathbb Z$ 上对应的差分矩阵，其是无穷阶的（
式 23 表明对应指标(j,m) 的矩阵元为 $\frac{(-1)^{j-m}}{(j-m)h}$，容易验证其是
Toeplitz 矩阵）：
\begin{equation}
h^{-1}
\begin{pmatrix}
&&&&&&&&\vdots&&&&&&&&\\
&&&&&&&&1/3&&&&&&&&\\
&&&&\ddots&&&&-1/2&&&&&&&&\\
&&&&\ddots&&&&1&&&&&&&&\\
&&&&\ddots&&&&0&&&&\ddots&&&&\\
&&&&&&&&-1&&&&\ddots&&&&\\
&&&&&&&&1/2&&&&\ddots&&&&\\
&&&&&&&&-1/3&&&&&&&&\\
&&&&&&&&\vdots&&&&&&
\end{pmatrix}.~
\end{equation}
周期情形
同样的，式 11 可以写为
\begin{equation}
\begin{aligned}
p(x)&=\frac{1}{2L}{\sum'}\limits_{k=-\frac{N}{2}}^{\frac{N}{2}}e^{ \mathrm{i} k\frac{\pi}{L}x}\hat u(k)\\
&=\frac{1}{2L}\sum_{m=1}^{N}u(x_m){\sum'}\limits_{k=-\frac{N}{2}}^{\frac{N}{2}}e^{ \mathrm{i} k\frac{\pi}{L}x}\hat \delta_m(k),~
\end{aligned}~
\end{equation}
而
\begin{equation}
\hat \delta_m(k)=h\sum_{i=1}^{N}e^{- \mathrm{i} k\frac{\pi}{L}x_i}\delta(x_i-x_m)=he^{- \mathrm{i} k\frac{\pi}{L}x_m}.~
\end{equation}
所以
\begin{equation}
\begin{aligned}
p(x)
&=\frac{h}{2L}\sum_{m=1}^{N}u(x_m){\sum'}\limits_{k=-\frac{N}{2}}^{\frac{N}{2}}e^{ \mathrm{i} k\frac{\pi}{L}(x-x_m)}\\
&=\frac{h}{2L}\sum_{m=1}^{N}u(x_m) \left(\frac{1}{2}\sum_{k=-\frac{N}{2}}^{\frac{N}{2}-1}e^{ \mathrm{i} k\frac{\pi}{L}(x-x_m)}+\frac{1}{2}\sum_{k=-\frac{N}{2}+1}^{\frac{N}{2}}e^{ \mathrm{i} k\frac{\pi}{L}(x-x_m)} \right) \\
&=\frac{h}{2L}\sum_{m=1}^{N}u(x_m)\cos \left(\frac{\pi}{2L}(x-x_m) \right) \sum_{k=-\frac{N}{2}+1/2}^{\frac{N}{2}-1/2}e^{ \mathrm{i} k\frac{\pi}{L}(x-x_m)}\\
&=\frac{h}{2L}\sum_{m=1}^{N}u(x_m)\cos \left(\frac{\pi}{2L}(x-x_m) \right) \frac{e^{ \mathrm{i} (-N/2+1/2)\frac{\pi}{L}(x-x_m)}-e^{ \mathrm{i} (N/2+1/2)\frac{\pi}{L}(x-x_m)}}{1-e^{ \mathrm{i} \frac{\pi}{L}(x-x_m)}}\\
&=\frac{h}{2L}\sum_{m=1}^{N}u(x_m)\cos \left(\frac{\pi}{2L}(x-x_m) \right) \frac{e^{- \mathrm{i} \frac{N\pi}{2L}(x-x_m)}-e^{ \mathrm{i} \frac{N\pi}{2L}(x-x_m)}}{e^{- \mathrm{i} \frac{\pi}{2L}(x-x_m)}-e^{ \mathrm{i} \frac{\pi}{2L}(x-x_m)}}\\
&=\frac{h}{2L}\sum_{m=1}^{N}u(x_m)\cos \left(\frac{\pi}{2L}(x-x_m) \right) \frac{\sin \left(\frac{N\pi}{2L}(x-x_m) \right) }{\sin \left(\frac{\pi}{2L}(x-x_m) \right) }\\
&=\sum_{m=1}^{N}u(x_m)\frac{\sin \left(\frac{N\pi}{2L}(x-x_m) \right) }{\frac{2L}{h}\tan \left(\frac{\pi}{2L}(x-x_m) \right) }.
\end{aligned}~
\end{equation}
由
\begin{equation}
\frac{\mathrm{d}{}}{\mathrm{d}{x}} \left(\frac{\sin \left(\frac{N\pi}{2L}(x-x_m) \right) }{\frac{2L}{h}\tan \left(\frac{\pi}{2L}(x-x_m) \right) } \right) |_{x=x_j}=\left\{\begin{aligned}
&0,\quad j=m\\
&\frac{\pi}{2L}(-1)^{(j-m)}\cot[\pi(j-m)/N],j\neq m
\end{aligned}\right.~
\end{equation}
从而
\begin{equation}
\begin{aligned}
w_j=p'(x_j)=\frac{\pi}{2L}\sum_{m=1,m\neq j}^{N}(-1)^{(j-m)}\cot[\pi(j-m)/N]u(x_m)
\end{aligned}~
\end{equation}
因此，差分矩阵为
\begin{equation}
\frac{\pi}{2L}
\begin{pmatrix}
0& \cot\left(\pi/N\right) &- \cot\left(2\pi/N\right) &\cdots& \cot\left((N-1)\pi/N\right) \\
- \cot\left(\pi/N\right) &\ddots&\ddots&\ddots&\ddots\\
\cot\left(2\pi/N\right) &\ddots&\ddots&\ddots&\ddots\\
\vdots&\ddots&\ddots&\ddots&\ddots\\
- \cot\left((N-1)\pi/N\right) &\ddots&\ddots&\ddots&\ddots\\
\end{pmatrix}.~
\end{equation}
致读者： 小时百科一直以来坚持所有内容免费无广告，这导致我们处于严重的亏损状态。 长此以往很可能会最终导致我们不得不选择大量广告以及内容付费等。 因此，我们请求广大读者
热心打赏 ，使网站得以健康发展。 如果看到这条信息的每位读者能慷慨打赏 20 元，我们一周就能脱离亏损， 并在接下来的一年里向所有读者继续免费提供优质内容。 但遗憾的是只有不到 1% 的读者愿意捐款， 他们的付出帮助了 99% 的读者免费获取知识， 我们在此表示感谢。