5.7 幂级数、Fourier 分析与正交展开

上一节讨论函数项级数的一般理论，本节处理两类最重要的函数项级数。第一类是幂级数 $\sum a_k{x}^k$，它是 Taylor 展开的严格版本：给定系数，级数在某个区间内定义出一个函数，这个函数自动具有任意阶导数，并且在收敛区间内可以逐项求导与逐项积分。第二类是三角级数，即 Fourier 级数，它把函数按正弦与余弦展开。两类展开有一个共同的框架：把函数看作无穷维向量，在内积空间里找一组正交基，展开系数就是函数在这组基上的投影。这个框架一旦建立，两类展开的平行结构就清楚了。

正交展开的实际价值来自计算。Fourier 级数的系数可以用快速 Fourier 变换（FFT）在 $O(N\log N)$ 的时间内算出来，这个算法是信号处理、数值求解偏微分方程与大规模数据分析的基础部件，本节实测朴素 DFT 与 FFT 的耗时差距（$N=4096$ 时相差四个数量级）。正交多项式（Legendre、Chebyshev）提供另一类基，它们的展开系数衰减速率直接反映函数的光滑程度，本节用三组数据对比几何收敛与代数收敛的差别。这两条线索在 5.6 节的逼近理论与 5.10 节的 $L^2$ 空间之间架了桥。

5.7.1 幂级数与解析函数

收敛半径

形如 $\sum _{k=0}^{\mathrm{\infty }}{a}_k(x-x_0{)}^k$ 的级数称幂级数（$x_0=0$ 时简写为 $\sum a_k{x}^k$）。幂级数的收敛集合总是以展开点为中心的一个区间，由 Cauchy-Hadamard 公式给出半径：

$$R=\frac{1}{\underset{k\to \mathrm{\infty }}{lim sup}|a_k{|}^{1/k}}$$

级数在 $|x-x_0|<R$ 内绝对收敛，在 $|x-x_0|>R$ 发散，在端点 $|x-x_0|=R$ 处的行为需要单独判断。比例形式（比值法）在极限存在时更便于计算：

$$R=\underset{k\to \mathrm{\infty }}{lim}\left|\frac{a_k}{a_{k+1}}\right|$$

Abel 定理说明收敛集合的内部性质：幂级数在收敛区间内的任意闭子区间 $[-r,r]$（$r<R$）上一致收敛（由 M 判别法，$|a_k{r}^k|$ 给出收敛的正项级数），因此和函数在 $(-R,R)$ 内连续。

逐项运算

在收敛区间内部，幂级数可以像多项式一样运算。

逐项求导：$f(x)=\sum a_k{x}^k$ 在 $|x|<R$ 内可导，且

$$f^{\prime }(x)=\sum _{k=1}^{\mathrm{\infty }}k{a}_k{x}^{k-1}$$

求导后的级数半径不变（$lim sup|k{a}_k{|}^{1/k}=lim sup|a_k{|}^{1/k}$，因为 $k^{1/k}\to 1$）。反复求导得 $a_k=\frac{f^{(k)}(0)}{k!}$，幂级数的系数由和函数在展开点的各阶导数唯一决定。逐项积分：${\int }_0^{x}f(t)dt=\sum \frac{a_k}{k+1}{x}^{k+1}$，半径同样不变。

这两条性质给出了求幂级数和函数的标准方法：把级数认成某个已知展开的导数或积分。例：

$$\sum _{k=1}^{\mathrm{\infty }}k{x}^k=x\sum _{k=1}^{\mathrm{\infty }}k{x}^{k-1}=x\cdot \frac{d}{dx}\left(\frac{1}{1-x}\right)=\frac{x}{(1-x{)}^2},|x|<1$$

例：$\sum _{k=1}^{\mathrm{\infty }}\frac{x^k}{k}=-\ln (1-x)$（对 $\sum x^k$ 逐项积分，$|x|<1$）。

逐项运算的条件

逐项求导需要 5.6 节的求导交换定理：导函数列在收敛区间内部一致收敛即满足。幂级数在内部的任意闭子区间上一致收敛，导数级数同样如此（半径相同），条件自动满足。这是幂级数比其他函数项级数好用的根本原因：把收敛半径算出来，逐项运算就都合法。端点处的行为需要单独讨论，端点上的逐项运算可能失效（$\sum \frac{x^k}{k}$ 在 $x=1$ 处发散，在 $x=-1$ 处条件收敛但不可逐项求导）。

Taylor 级数与解析性

$f$ 在 $x_0$ 处的 Taylor 级数是 $\sum \frac{f^{(k)}(x_0)}{k!}(x-x_0{)}^k$。5.4 节的 Taylor 公式说明有限项的 Taylor 多项式是好的局部近似，Taylor 级数是这个多项式在 $n\to \mathrm{\infty }$ 时的极限。级数是否收敛到 $f$ 本身依赖余项的行为：$f$ 在 $x_0$ 附近实解析指存在 $\delta >0$ 使得 $|x-x_0|<\delta$ 时 Taylor 级数收敛到 $f(x)$。

解析与任意阶可导不等价。标准例子是

$$f(x)=\{\begin{array}{ll}{e}^{-1/x^2},& x\ne 0\\ 0,& x=0\end{array}$$

这个函数在 $\mathbb{R}$ 上任意阶可导，且所有导数在 $x=0$ 处为零（$e^{-1/x^2}$ 在零点附近趋于零的速度快于任何幂函数），因此 Taylor 级数是零级数，只在 $x=0$ 处与 $f$ 相等。这个例子说明「有 Taylor 级数」与「等于 Taylor 级数」是两件事。复变函数的情形不同：复平面上一次可导蕴含任意阶可导，也蕴含解析，因此复解析函数的 Taylor 级数一定收敛到自身，收敛半径等于展开点到最近奇点的距离。Runge 函数 $f(x)=\frac{1}{1+25x^2}$ 在实轴上处处解析，但 Taylor 级数在 $|x|>\frac{1}{5}$ 处发散，原因就是复平面上的奇点 $x=\pm \frac{i}{5}$ 限制了半径，这一点在 5.6 节的插值讨论中已经出现过。

实验一用三组系数验证 Cauchy-Hadamard 公式。$a_k=\frac{1}{k!}$ 时 $|a_k{|}^{1/k}$ 与比值都趋于零，半径无穷；$a_k=\frac{1}{2^k}$ 时两者恒为 $\frac{1}{2}$，半径 $2$；$a_k=\frac{k^k}{k!}$ 时两个量都趋于欧拉数 $e=2.71828$，半径 $\frac{1}{e}=0.367879$。第三组的数值收敛较慢（$k=100$ 时比值是 $2.704814$，还差 $0.5\mathrm{\%}$），因为它以 ${(1+\frac{1}{k})}^k$ 的速率逼近 $e$。

5.7.2 内积空间与正交展开

从向量到函数

第一章的内积把「夹角」与「投影」引入有限维向量空间。函数之间也可以定义内积，取 $L^2$ 内积（以区间 $[a,b]$ 为底）：

$$⟨f,g⟩={\int }_a^{b}f(x)g(x)dx,\parallel f{\parallel }_2=\sqrt{⟨f,f⟩}$$

对应的范数 $\parallel \cdot {\parallel }_2$ 称 $L^2$ 范数，$\parallel f-g{\parallel }_2$ 是两条曲线之间的均方距离。这一点与 5.6 节的上确界范数不同：上确界范数度量最大竖直距离，均方距离度量平均意义下的偏差。两套范数给出两种收敛，一致收敛与均方收敛，后者更弱但更容易达到，Fourier 级数的收敛性在均方意义下最干净。

若 $⟨{\varphi }_i,{\varphi }_j⟩=0$（$i\ne j$），称函数系 $\{{\varphi }_i\}$ 正交；若还有 $\parallel {\varphi }_i{\parallel }_2=1$，称标准正交。正交系的核心好处是系数的解耦：要求 $\sum c_k{\varphi }_k$ 逼近 $f$，最优系数可以直接写成投影

$$c_k=\frac{⟨f,{\varphi }_k⟩}{⟨{\varphi }_k,{\varphi }_k⟩}$$

不需要解线性方程组。原因是一元二次函数 $\parallel f-c{\varphi }_k{\parallel }_2^2$ 对 $c$ 求导时其他项都消失（正交性），方程退化为一次方程。

最佳均方逼近与 Bessel 不等式

设 $\{{\varphi }_k\}$ 是正交系，$S_n=\sum _{k=1}^n{c}_k{\varphi }_k$ 是用前 $n$ 个基函数给出的部分和。在 $n$ 维子空间 $\mathrm{span}\{{\varphi }_1,\dots ,{\varphi }_n\}$ 中，$S_n$ 是 $f$ 的最佳均方逼近：任何其他线性组合的均方误差都不小于 $\parallel f-S_n{\parallel }_2$。这条结论把逼近问题化归为计算投影，也是最小二乘法的几何内容（第一章的正交投影一节）。

由 $\parallel f-S_n{\parallel }_2^2=\parallel f{\parallel }_2^2-\sum _{k=1}^n|c_k{|}^2\parallel {\varphi }_k{\parallel }_2^2\ge 0$ 得 Bessel 不等式

$$\sum _{k=1}^{\mathrm{\infty }}|c_k{|}^2\parallel {\varphi }_k{\parallel }_2^2\le \parallel f{\parallel }_2^2$$

左边是展开系数的能量，右边是函数本身的能量。推论：$c_k\to 0$，系数必须衰减。这条推论给出了必要条件的定量形式，但没有说明衰减多快，衰减速率依赖函数的光滑程度（5.7.5 节）。

完备性与 Parseval 等式

若对每个 $f$ 都有 $\parallel f-S_n{\parallel }_2\to 0$，称正交系完备。完备时 Bessel 不等式取等号，得到 Parseval 等式

$$\sum _{k=1}^{\mathrm{\infty }}|c_k{|}^2\parallel {\varphi }_k{\parallel }_2^2=\parallel f{\parallel }_2^2$$

Parseval 等式的读法是能量守恒：函数在原始空间的能量等于它在各个正交方向上的投影能量之和。这条等式把函数的积分（连续量）化成级数求和（离散量），在概率论中对应方差分解（第二章的方差与条件期望），在信号处理中对应功率谱。三角多项式系与多项式系（Legendre、Chebyshev 等）都是完备的正交系，各自给出一个具体的 Parseval 等式。

5.7.3 Fourier 级数

系数与基本性质

在 $[-\pi ,\pi ]$ 上，三角函数系

$$\{\frac{1}{\sqrt{2\pi }},\frac{\cos x}{\sqrt{\pi }},\frac{\sin x}{\sqrt{\pi }},\dots ,\frac{\cos kx}{\sqrt{\pi }},\frac{\sin kx}{\sqrt{\pi }},\dots \}$$

是标准正交系（由三角函数的正交关系 ${\int }_{-\pi }^{\pi }\cos kx\cos lxdx=\pi {\delta }_{kl}$ 等得到）。函数 $f$ 的 Fourier 级数是

$$f(x)\sim \frac{a_0}{2}+\sum _{k=1}^{\mathrm{\infty }}(a_k\cos kx+b_k\sin kx)$$

系数由投影给出：

$$a_k=\frac{1}{\pi }{\int }_{-\pi }^{\pi }f(x)\cos kxdx,b_k=\frac{1}{\pi }{\int }_{-\pi }^{\pi }f(x)\sin kxdx$$

三条基本性质决定了展开的形态。奇偶性简化：偶函数的 $b_k$ 全为零（只剩余弦项），奇函数的 $a_k$ 全为零（只剩正弦项）。频率分解：Fourier 展开把一个函数分解成不同频率的振荡之和，$a_k$、$b_k$ 构成频谱，这是它在信号处理与谱分析中不可替代的原因。平移规律：把 $f$ 平移，系数的绝对值不变而相位按频率线性改变，这是卷积定理与滤波器设计的基础。

收敛性

Fourier 级数的收敛有三种含义，强弱不同。

均方收敛：对任何 $f\in L^2[-\pi ,\pi ]$，部分和的均方误差 $\parallel f-S_n{\parallel }_2\to 0$，Parseval 等式成立。这条结论对平方可积函数都成立，条件最弱。

逐点收敛与一致收敛：若 $f$ 在 $[-\pi ,\pi ]$ 上除有限个第一类间断点外连续、在其余点可导，则级数在每个连续点收敛到 $f(x)$，在间断点 $x_0$ 处收敛到 $\frac{f(x_0^{+})+f(x_0^{-})}{2}$（左右极限的算术平均）。若 $f$ 连续、分段 $C^1$ 且满足 $f(-\pi )=f(\pi )$（周期延拓后连续），则级数一致收敛。收敛速率由光滑程度决定：第 $k$ 项系数的量级是 $O(k^{-p})$，其中 $p$ 与函数的可导阶数相关（函数 $p-1$ 阶导数有界则系数 $O(k^{-p})$）。方波的系数是 $O(1/k)$（只有间断），三角波的系数是 $O(1/k^2)$（连续但导数有跳跃），两者的收敛性质差别很大。

发散的存在：存在连续函数的 Fourier 级数在某点发散（du Bois-Reymond 构造），因此「连续」不足以保证逐点收敛。均方收敛则没有这个问题，这也是 $L^2$ 理论成为 Fourier 分析主流框架的原因。

Gibbs 现象

方波 $f(x)=\mathrm{sgn}(x)$ 的 Fourier 级数是

$$f(x)\sim \frac{4}{\pi }\sum _{k=1,3,5,\dots }^{\mathrm{\infty }}\frac{\sin kx}{k}$$

实验二给出了部分和的过冲数据。$N=5$ 时最大过冲 $1.188357$，$N=501$ 时 $1.178980$，过冲比例稳定在 $17.9\mathrm{\%}$ 附近（相对跳跃幅度 $2$ 的比例是 $8.949\mathrm{\%}$）。过冲不随 $N$ 消失：极限值是

$$\underset{N\to \mathrm{\infty }}{lim}max{S}_N=\frac{2}{\pi }{\int }_0^{\pi }\frac{\sin t}{t}dt=1.178980$$

即 $2/\pi \cdot \mathrm{Si}(\pi )$，其中 $\mathrm{Si}$ 是正弦积分。这个值比 $1$ 高出 $0.178980$，与 $N$ 无关。过冲的位置随 $N$ 移动：极值点在 $x\approx \frac{\pi }{N}$ 附近，因此过冲的宽度趋于零。这解释了 Gibbs 现象的完整图景：函数值在几乎处处收敛，但间断点两侧出现宽度趋零、幅度固定的振荡，最大偏差不趋于零，一致收敛不成立（这与 $x^n$ 的情形在结构上属于同一类）。Gibbs 现象出现在图像压缩、滤波器设计中就形成振铃（ringing）伪影，处理方法包括加窗（用平滑因子压低高频系数）或用总变差型正则化。

作为对照，三角波（连续但导数有跳跃）的系数按 $\frac{1}{k^2}$ 衰减，部分和的最大值从 $1.414711$ 单调升向 $\frac{\pi }{2}=1.570796$，不出现过冲，一致收敛成立。两段数据的对照说明：间断造成过冲，连续性与更高的可导阶数带来更好的收敛。

均方误差与收敛速率

方波的部分和均方误差从 $N=5$ 的 $0.648568$ 降到 $N=501$ 的 $0.071333$，按 $\frac{1}{\sqrt{N}}$ 衰减。这个速率可以直接算出：由 Parseval 等式，尾部的能量是

$$\parallel f-S_N{\parallel }_2^2=\pi \sum _{k>N,k\text{奇}}{\left(\frac{4}{\pi k}\right)}^2=\frac{16}{\pi }\sum _{k>N,k\text{奇}}\frac{1}{k^2}\approx \frac{8}{\pi N}$$

因此 $\parallel f-S_N{\parallel }_2\approx \sqrt{\frac{8}{\pi N}}$。代入 $N=501$ 得 $0.0713$，与实测的 $0.071333$ 一致；$N=101$ 得 $0.1588$，与实测的 $0.158052$ 一致。三角波的对照数据：$N=51$ 时均方误差 $0.002457$，因为它按 $N^{-3/2}$ 衰减（系数按 $k^{-2}$，尾部能量的系数按 $k^{-4}$，求和得 $N^{-3}$，开方得 $N^{-3/2}$）。两条速率公式说明：均方误差也是由系数的衰减速率决定的，而系数衰减由光滑程度决定。

三种收敛的关系

在有限区间上，一致收敛蕴含均方收敛（$\parallel f{\parallel }_2\le \sqrt{b-a}\parallel f{\parallel }_{\mathrm{\infty }}$），但均方收敛不蕴含逐点收敛（改变函数在零测集上的值不影响均方误差，却可能破坏每个点上的收敛）。方波的部分和在均方意义下收敛，在间断点处按算术平均收敛，在 Gibbs 过冲处的取值不收敛到函数值。三套收敛各自回答不同的问题，讨论收敛性时应当先说明所用的是哪一种。

5.7.4 离散 Fourier 变换与快速算法

DFT 的定义与性质

实际数据是有限个采样点，连续 Fourier 变换的对应物是离散 Fourier 变换（DFT）。对长度 $N$ 的序列 $x_0,\dots ,x_{N-1}$（可以视为复向量），定义

$$X_k=\sum _{n=0}^{N-1}{x}_n{e}^{-2\pi ikn/N},k=0,1,\dots ,N-1$$

$X_k$ 是复数，$|X_k|$ 是频率 $k$ 的分量幅度，$\arg X_k$ 是相位。逆变换由

$$x_n=\frac{1}{N}\sum _{k=0}^{N-1}{X}_k{e}^{2\pi ikn/N}$$

给出（正逆变换的差别只有符号与系数 $\frac{1}{N}$）。DFT 的矩阵形式是 $X=Wx$，其中 $W$ 的元素是 $W_{kn}=e^{-2\pi ikn/N}$，$W$ 是酉矩阵（$W^{\ast }W=NI$），这一条性质给出离散版本的 Parseval 等式

$$\sum _{k=0}^{N-1}|X_k{|}^2=N\sum _{n=0}^{N-1}|x_n{|}^2$$

三条常用性质使 DFT 成为计算工具。卷积定理：序列的循环卷积对应频谱的逐点相乘，因此滤波运算可以在频域完成（代价从 $O(N^2)$ 降到 $O(N\log N)$）。平移与调制：时域的平移对应频域的相位因子，频域的平移对应时域的调制。实信号的共轭对称：$x_n$ 为实数时 $X_{N-k}=\overline{X_k}$，频谱的一半信息是冗余的，存储与计算可以只保留一半。实验三的最后一个例子展示了这条性质：信号 $1.5\sin 3t+0.8\cos 7t+0.3$ 的频谱在 $k=3$ 与 $k=61$（即 $N-3$）处成对出现，幅度都是 $0.750000$；$k=7$ 与 $k=57$ 成对出现，幅度 $0.400000$；直流分量在 $k=0$ 处为 $0.300000$。

快速 Fourier 变换

DFT 的直接计算需要 $N^2$ 次复数乘法，$N={10}^6$ 时约 ${10}^{12}$ 次运算，不可接受。快速 Fourier 变换（FFT）利用 $e^{-2\pi ikn/N}$ 的周期性把计算量降到 $O(N\log N)$，方法是分治：把长度为 $N$ 的变换拆成两个长度为 $\frac{N}{2}$ 的变换（按偶数下标与奇数下标分组），

$$X_k=E_k+e^{-2\pi ik/N}{O}_k,X_{k+N/2}=E_k-e^{-2\pi ik/N}{O}_k$$

其中 $E$、$O$ 分别是偶下标与奇下标的子序列的 DFT。递归下去，每层的工作量是 $O(N)$，共 ${\log }_{2}N$ 层，总复杂度 $O(N\log N)$。最高效的实现（Cooley-Tukey 算法）对 $N$ 为 $2$ 的幂的情形递归拆分，并且用原地置换与蝶形运算减少存储访问；一般 $N$ 的实现先用 Bluestein 或 Rader 算法把问题化归为素数长度的变换。

实验三的耗时数据把差距量化了（时间依赖机器，趋势是稳定的）：朴素 DFT 的耗时随 $N$ 每翻一倍增加到 $4$ 倍（$N=256$ 时 $0.001844$ 秒，$N=4096$ 时 $0.406679$ 秒，$16$ 倍数据量带来 $220$ 倍耗时，接近 $N^2$ 的 $256$ 倍），FFT 的耗时增长缓慢（$N=4096$ 时 $0.000029$ 秒），两者的比值从 $283.7$ 扩大到 $14023.4$。两者的结果一致（最大偏差 ${10}^{-12}$ 到 ${10}^{-10}$ 量级，属于浮点舍入），说明 FFT 只是重新组织了运算次序，没有改变数学内容。

FFT 的实用地位来自它的复杂度：$N={10}^6$ 时 $O(N\log N)\approx 2\times {10}^7$ 次运算，毫秒级完成；卷积定理把滤波、相关、多项式乘法都化为 FFT 加逐点乘法，这条路线在音频处理、图像处理、雷达与通信、数值偏微分方程（谱方法）中都是标准做法。

5.7.5 正交多项式

Legendre 与 Chebyshev 多项式

在 $[-1,1]$ 上取权函数 $w(x)\equiv 1$，对单项式 $1,x,x^2,\dots$ 做 Gram-Schmidt 正交化，得到的正交多项式系是 Legendre 多项式 $\{P_k\}$：

$$P_0=1,P_1=x,P_2=\frac{3x^2-1}{2},P_3=\frac{5x^3-3x}{2},\dots$$

满足 ${\int }_{-1}^1{P}_i(x)P_j(x)dx=\frac{2}{2i+1}{\delta }_{ij}$。Legendre 多项式的零点正是 5.5 节 Gauss-Legendre 求积的节点，那里提到的「节点与权重需要预计算」在这里找到了数学来源：最优节点是正交多项式的零点。

取权函数 $w(x)=\frac{1}{\sqrt{1-x^2}}$，得到的正交系是 Chebyshev 多项式 $\{T_k\}$，定义为 $T_k(\cos \theta )=\cos k\theta$，满足

$${\int }_{-1}^1\frac{T_i(x)T_j(x)}{\sqrt{1-x^2}}dx=\{\begin{array}{ll}\pi ,& i=j=0\\ \frac{\pi }{2},& i=j>0\\ 0,& i\ne j\end{array}$$

$T_k$ 的零点与极值点正是 5.6 节插值讨论中的 Chebyshev 节点，$\underset{[-1,1]}{max}|T_k|=1$ 这条性质给出了节点多项式的最小最大界 $2^{1-n}$。两条线索在这里接上：插值节点、正交多项式的零点、Gauss 求积的节点，是同一个数学对象在不同问题中的三个侧面。

展开系数的衰减与光滑程度

把函数按正交多项式展开，$f\sim \sum c_k{\varphi }_k$，$c_k=\frac{⟨f,{\varphi }_k⟩}{⟨{\varphi }_k,{\varphi }_k⟩}$。系数的衰减速率由 $f$ 的光滑程度决定，这一条结论是正交展开最重要的实用内容：

解析函数：系数按几何速率衰减，$|c_k|\le C{\rho }^{-k}$（$\rho$ 由奇点到区间的距离决定，与 5.6 节的 Bernstein 椭圆参数相同），截断误差也按几何速率下降。

$p$ 阶连续可导的函数：系数按 $|c_k|\le \frac{C}{k^p}$ 衰减，截断误差按 $n^{-(p-1)}$ 量级下降。

只有连续性（含拐角）：系数按 $O(k^{-2})$ 量级衰减，截断误差按 $O(n^{-1})$ 下降。

实验四用三组数据验证这个阶梯。$f(x)=e^x$ 在整个复平面解析，Chebyshev 部分和的误差从 $n=2$ 的 $4.99\times {10}^{-2}$ 降到 $n=16$ 的 $1.78\times {10}^{-15}$（机器精度），之后不再下降（$n=64$ 时是 $1.60\times {10}^{-14}$，舍入噪声主导），四次加倍共降低十三个数量级。$f(x)=\left|x\right|$ 在 $x=0$ 处有拐角，误差从 $2.12\times {10}^{-1}$ 降到 $9.80\times {10}^{-3}$，$n$ 每翻一倍误差减半，对应 $O(n^{-1})$。$f(x)=\frac{1}{1+25x^2}$ 解析但有距离 $0.2$ 的复奇点，误差从 $5.40\times {10}^{-1}$ 降到 $2.41\times {10}^{-6}$，比值 $0.034$ 对应 ${\rho }^{-2}$（$\rho =1.2198$，${\rho }^{-2}=0.672$；实测比值比理论值小，因为误差在不同位置取到最大值，几何收敛的渐近阶段在 $n$ 较大时才充分建立），$n=64$ 时的理论量级 ${\rho }^{-64}={1.22}^{-64}\approx 3.0\times {10}^{-6}$ 与实测的 $2.41\times {10}^{-6}$ 相符。

Legendre 展开（用均匀网格上的最小二乘拟合近似）给出同样的阶梯：$e^x$ 在次数 $16$ 达到机器精度（$5.83\times {10}^{-14}$），$\left|x\right|$ 按 $O(n^{-1})$ 衰减（$1.17\times {10}^{-1}\to 1.90\times {10}^{-2}$），$\frac{1}{1+25x^2}$ 按几何速率下降（$3.30\times {10}^{-1}\to 1.36\times {10}^{-3}$）。两套基的差别在常数因子：对同一函数，Chebyshev 展开的最大误差通常略小于 Legendre 展开（因为 Chebyshev 基的误差分布更均匀），而 Legendre 基的展开系数与 Gauss-Legendre 求积直接对应，在数值积分中更方便。

实用结论：看到展开系数按代数速率衰减，说明函数有间断或拐角，提高截断阶数收益有限（要降低一个数量级的误差，需要提高一个数量级的阶数）；看到几何衰减，说明函数解析，提高阶数的收益是指数级的。这条判据在谱方法、压缩感知与模型逼近中被用作「选择基函数与截断阶数」的依据。

5.7.6 数值实验

实验一：收敛半径的经验测定

对三组系数计算 $|a_k{|}^{1/k}$ 与 $\left|\frac{a_{k+1}}{a_k}\right|$ 两个估计量，对照理论半径。

import math

def ratios(a, k):
    """返回 (|a_k|^{1/k}, |a_{k+1}/a_k|)"""
    r1 = abs(a(k)) ** (1.0 / k)
    r2 = abs(a(k + 1) / a(k))
    return r1, r2

series = [
    ("a_k = 1/k!", lambda k: 1.0 / math.factorial(k), "理论半径 无穷大"),
    ("a_k = 1/2^k", lambda k: 1.0 / 2.0 ** k, "理论半径 2"),
    ("a_k = k^k/k!", lambda k: k ** k / math.factorial(k), "理论半径 1/e = 0.367879"),
]
for name, a, theory in series:
    print(f"\n{name}    （{theory}）")
    print(f"{'k':>6} {'|a_k|^(1/k)':>16} {'|a_(k+1)/a_k|':>18}")
    for k in [5, 10, 20, 50, 100]:
        r1, r2 = ratios(a, k)
        print(f"{k:>6} {r1:>16.6f} {r2:>18.6f}")

print("部分和检验（递推计算，避免阶乘溢出）")
print(f"{'x':>8} {'N=10':>16} {'N=50':>16} {'N=200':>16} {'e^(-x) 精确值':>16}")
for x in [-2.0, -1.0, 0.5, 1.0, 2.0]:
    vals = []
    for N in [10, 50, 200]:
        term, s = 1.0, 1.0
        for k in range(1, N + 1):
            term *= (-x) / k
            s += term
        vals.append(s)
    print(f"{x:>8.1f} {vals[0]:>16.10f} {vals[1]:>16.10f} {vals[2]:>16.10f} {math.exp(-x):>16.10f}")

print("几何级数 sum x^k 在半径 1 两侧")
for x in [0.5, 0.9, 0.99, 1.0, 1.01, 1.1]:
    s = sum(x ** k for k in range(101))
    exact = 1.0 / (1.0 - x) if abs(x - 1) > 1e-15 else float("inf")
    print(f"  x = {x:>5.2f}: 部分和(N=100) = {s:>20.8f}   1/(1-x) = {exact:>20.8f}")

输出（节选主要部分）：

a_k = 1/k!    （理论半径 无穷大）
     k      |a_k|^(1/k)      |a_(k+1)/a_k|
     5         0.383852           0.166667
    10         0.220813           0.090909
    20         0.120419           0.047619
    50         0.051326           0.019608
   100         0.026321           0.009901

a_k = 1/2^k    （理论半径 2）
     k      |a_k|^(1/k)      |a_(k+1)/a_k|
     5         0.500000           0.500000
    10         0.500000           0.500000
    20         0.500000           0.500000
    50         0.500000           0.500000
   100         0.500000           0.500000

a_k = k^k/k!    （理论半径 1/e = 0.367879）
     k      |a_k|^(1/k)      |a_(k+1)/a_k|
     5         1.919260           2.488320
    10         2.208125           2.593742
    20         2.408373           2.653298
    50         2.566306           2.691588
   100         2.632085           2.704814

部分和检验（递推计算，避免阶乘溢出）
       x             N=10             N=50            N=200       e^(-x) 精确值
    -2.0     7.3889947090     7.3890560989     7.3890560989     7.3890560989
    -1.0     2.7182818011     2.7182818285     2.7182818285     2.7182818285
     0.5     0.6065306597     0.6065306597     0.6065306597     0.6065306597
     1.0     0.3678794643     0.3678794412     0.3678794412     0.3678794412
     2.0     0.1353791887     0.1353352832     0.1353352832     0.1353352832

几何级数 sum x^k 在半径 1 两侧
  x =  0.50: 部分和(N=100) =           2.00000000   1/(1-x) =           2.00000000
  x =  0.90: 部分和(N=100) =           9.99976095   1/(1-x) =          10.00000000
  x =  0.99: 部分和(N=100) =          63.76279821   1/(1-x) =         100.00000000
  x =  1.00: 部分和(N=100) =         101.00000000   1/(1-x) =                  inf
  x =  1.01: 部分和(N=100) =         173.18619677   1/(1-x) =        -100.00000000
  x =  1.10: 部分和(N=100) =      151576.73573805   1/(1-x) =         -10.00000000

三组系数给出三种半径的形态。第一组的两个估计量都单调趋于零（$|a_k{|}^{1/k}$ 从 $0.383852$ 降到 $0.026321$，比值从 $0.166667$ 降到 $0.009901$），倒数趋于无穷，半径无穷大。第二组的两个估计量恒为 $0.500000$ 与 $n$ 无关，半径 $2$，精确值。第三组的两个估计量都缓慢上升并逼近 $e=2.718282$：$|a_k{|}^{1/k}$ 从 $1.919260$ 升到 $2.632085$，比值从 $2.488320$ 升到 $2.704814$，$k=100$ 时还差 $0.5\mathrm{\%}$。逼近速度慢的原因可以在公式里看出：

$$\frac{a_{k+1}}{a_k}=\frac{(k+1{)}^{k+1}/(k+1)!}{k^k/k!}={(1+\frac{1}{k})}^k\to e$$

这个量以 $1-\frac{1}{2k}$ 的相对速率逼近 $e$（$k=100$ 时是 $1-0.005$，与数据中 $\frac{2.704814}{2.718282}=0.99505$ 相符）。半径 $R=\frac{1}{e}=0.367879$。

部分和检验一行数据说明「半径无穷」的实际含义：$x=\pm 2$ 时 $N=50$ 的部分和已经与 $e^{\mp 2}$ 的十位小数完全一致，$e^{-2}=0.1353352832$ 与 $0.1353352832$ 逐位相同。$N=10$ 时 $x=2$ 的部分和 $0.1353791887$ 与精确值的差是 $4.4\times {10}^{-5}$，这一项来自被截断的 $k=11$ 项（$\frac{2^{11}}{11!}=4.4\times {10}^{-5}$），与公式一致。

几何级数的六行数据把「半径内外」的差别演示得很直接。$x=0.5$ 时部分和与 $\frac{1}{1-x}$ 完全一致（${0.5}^{101}=4\times {10}^{-31}$，在双精度下不可见）；$x=0.9$ 时误差 $2.39\times {10}^{-4}$，等于尾部 $\frac{{0.9}^{101}}{0.1}$；$x=0.99$ 时部分和只有 $63.76$，离 $\frac{1}{0.01}=100$ 还差 $36\mathrm{\%}$，因为收敛速率是 $0.99$ 的幂，$N=100$ 只把误差降到 ${0.99}^{101}/(1-0.99)=36.2$；$x=1$ 时部分和等于 $101$，发散；$x=1.01$ 与 $x=1.1$ 时部分和是正的大数（$173.19$ 与 $151576.74$），而解析公式 $\frac{1}{1-x}$ 给出负数（$-100$ 与 $-10$）。后面两行说明解析公式在半径外不能使用：$\frac{1}{1-x}$ 的解析延拓在其他位置有定义，但与级数的和无关。

实验二：Fourier 级数与 Gibbs 现象

计算方波与三角波的 Fourier 部分和，观察最大过冲与均方误差的收敛。

import math
import numpy as np

xg = np.linspace(-math.pi + 1e-9, math.pi - 1e-9, 400001)
f = np.where(xg >= 0, 1.0, -1.0)
dx = xg[1] - xg[0]
print(f"{'N':>6} {'最大过冲':>12} {'过冲比例':>12} {'相对跳跃的过冲':>16} {'L2 误差':>13}")
for N in [5, 11, 51, 101, 501]:
    s = np.zeros_like(xg)
    for k in range(1, N + 1, 2):
        s += (4.0 / math.pi) * np.sin(k * xg) / k
    over = float(np.max(s))
    l2 = math.sqrt(float(np.sum((s - f) ** 2) * dx))
    print(f"{N:>6} {over:>12.6f} {(over - 1.0):>12.4%} {(over - 1.0) / 2.0:>16.4%} {l2:>13.6f}")

tri = np.abs(xg) - math.pi / 2
print(f"{'N':>6} {'max|S_N|':>12} {'L2 误差':>13}")
for N in [3, 7, 15, 51]:
    s = np.zeros_like(xg)
    for k in range(1, N + 1, 2):
        s += -(4.0 / math.pi) * np.cos(k * xg) / (k * k)
    l2 = math.sqrt(float(np.sum((s - tri) ** 2) * dx))
    print(f"{N:>6} {float(np.max(np.abs(s))):>12.6f} {l2:>13.6f}")

输出：

     N         最大过冲         过冲比例          相对跳跃的过冲         L2 误差
     5     1.188357     18.8357%          9.4178%      0.648568
    11     1.181302     18.1302%          9.0651%      0.460148
    51     1.179103     17.9103%          8.9552%      0.221315
   101     1.179012     17.9012%          8.9506%      0.158052
   501     1.178980     17.8980%          8.9490%      0.071333
Wilbraham-Gibbs 常数的理论极限：(2/pi) * int_0^pi sin(t)/t dt = 1.178980
理论最大过冲限值：0.178980
     N     max|S_N|         L2 误差
     3     1.414711      0.108991
     7     1.491625      0.040109
    15     1.531059      0.014341
    51     1.558555      0.002457

方波的五行数据把 Gibbs 现象的三个特征都写出来了。过冲的幅度趋于一个非零极限：$N=5$ 时 $1.188357$，$N=501$ 时 $1.178980$，极限值 $1.178980$ 与理论值 $\frac{2}{\pi }\mathrm{Si}(\pi )$ 逐位相符（$\mathrm{Si}(\pi )=1.851937051982467$，乘 $\frac{2}{\pi }$ 得 $1.178980$）。过冲不消失，这是与「收敛」直接冲突的一点：部分和的最大值不趋于 $1$，一致收敛不成立。过冲的相对幅度稳定在 $8.95\mathrm{\%}$（相对跳跃 $2$），或用 $17.9\mathrm{\%}$ 相对于峰值 $1$ 表示，两个数字是同一个量的两种写法。

均方误差的走向与过冲相反。它按预期收敛：$0.648568\to 0.460148\to 0.221315\to 0.158052\to 0.071333$。用理论公式 $\sqrt{\frac{8}{\pi N}}$ 核对：$N=101$ 时 $\sqrt{\frac{8}{3.14159\times 101}}=0.1588$（实测 $0.158052$），$N=501$ 时 $\sqrt{\frac{8}{3.14159\times 501}}=0.0713$（实测 $0.071333$）。两个量级上的偏差来自公式是渐近估计（忽略了 $N$ 以下奇数次项的精确求和与端点效应）。这组数据说明均方收敛与一致收敛的差别是实质性的：能量意义上的误差趋于零，最大偏差意义上的误差不趋于零。

三角波的四行数据是对照。$N$ 从 $3$ 增到 $51$，$max|S_N|$ 从 $1.414711$ 升到 $1.558555$，单调趋近 $\frac{\pi }{2}=1.570796$，从下方逼近而且没有越过极限（不存在过冲）。均方误差从 $0.108991$ 降到 $0.002457$，衰减比方波快得多（$N$ 增大约 $17$ 倍，误差降到 $\frac{1}{44}$，对应 $N^{-3/2}$ 的 ${17}^{1.5}=70$ 倍量级）。两条曲线的差别全部来自连续性与光滑程度：方波有跳跃（系数 $O(1/k)$），三角波连续但导数有跳跃（系数 $O(1/k^2)$）。

实验三：朴素 DFT 与 FFT 的耗时

实现 $O(N^2)$ 的朴素 DFT 并与 numpy.fft 对照，同时验证结果一致性与频谱的对称性。耗时取三次运行的最小值，减少系统负载带来的波动。

import math
import time
import numpy as np

def naive_dft(x):
    N = len(x)
    n = np.arange(N)
    k = n.reshape((N, 1))
    return np.exp(-2j * np.pi * k * n / N) @ x

def time_call(fn, *args):
    t0 = time.perf_counter()
    fn(*args)
    return time.perf_counter() - t0

rng = np.random.default_rng(0)
print(f"{'N':>6} {'朴素 DFT (s)':>16} {'FFT (s)':>14} {'比值':>10} {'最大偏差':>14}")
for N in [256, 512, 1024, 2048, 4096]:
    x = rng.normal(size=N)
    tn = min(time_call(naive_dft, x) for _ in range(3))
    tf = min(time_call(np.fft.fft, x) for _ in range(3))
    Xn = naive_dft(x)
    Xf = np.fft.fft(x)
    print(f"{N:>6} {tn:>16.6f} {tf:>14.6f} {tn / max(tf, 1e-9):>10.1f} "
          f"{float(np.max(np.abs(Xn - Xf))):>14.2e}")

N = 64
t = np.arange(N) * 2 * math.pi / N
sig = 1.5 * np.sin(3 * t) + 0.8 * np.cos(7 * t) + 0.3
X = np.fft.fft(sig) / N
print("非零频率分量（幅度大于 0.05 的）:")
for k in range(N):
    if abs(X[k]) > 0.05:
        freq = k if k <= N // 2 else k - N
        print(f"  k = {k:>3} (频率 {freq:>3}): 幅度 = {abs(X[k]):.6f}, "
              f"相位 = {math.atan2(X[k].imag, X[k].real):+.6f}")

输出：

     N       朴素 DFT (s)        FFT (s)         比值           最大偏差
   256         0.001844       0.000007      283.7       3.68e-12
   512         0.007328       0.000008      939.5       7.99e-12
  1024         0.027501       0.000010     2696.2       2.70e-11
  2048         0.107505       0.000015     6935.8       6.71e-11
  4096         0.406679       0.000029    14023.4       2.29e-10
非零频率分量（幅度大于 0.05 的）:
  k =   0 (频率   0): 幅度 = 0.300000, 相位 = +0.000000
  k =   3 (频率   3): 幅度 = 0.750000, 相位 = -1.570796
  k =   7 (频率   7): 幅度 = 0.400000, 相位 = -0.000000
  k =  57 (频率  -7): 幅度 = 0.400000, 相位 = +0.000000
  k =  61 (频率  -3): 幅度 = 0.750000, 相位 = +1.570796

耗时数据把复杂度差别量化了。朴素 DFT 的耗时随 $N$ 翻倍增至约 $4$ 倍：$0.001844\to 0.007328\to 0.027501\to 0.107505\to 0.406679$，相邻比值是 $3.97$、$3.75$、$3.91$、$3.78$，与 $N^2$ 的 $4$ 倍相符（偏离来自缓存与 BLAS 调用的常数开销）。FFT 的耗时从 $0.000007$ 秒增到 $0.000029$ 秒，$16$ 倍数据量只带来 $4$ 倍耗时，接近 $N\log N$ 的 $16\times \frac{12}{8}=24$ 倍（实际比这个估计更快，因为 $N$ 小时的部分开销被摊薄）。两者的比值从 $283.7$ 扩大到 $14023.4$，差距随 $N$ 持续扩大。结果的最大偏差在 ${10}^{-12}$ 到 ${10}^{-10}$ 之间，与 FFT 与朴素 DFT 所用求和次序不同带来的舍入误差量级一致（$N=4096$ 时约 $2.3\times {10}^{-10}$，相对量级 ${10}^{-13}$），说明 FFT 在数值上与朴素 DFT 等价。

频谱的例子验证了实信号 DFT 的共轭对称与幅度公式。输入 $1.5\sin 3t+0.8\cos 7t+0.3$ 的解析频谱是：$\sin$ 项给出成对的共轭系数，频率 $3$ 处幅度 $\frac{1.5}{2}=0.75$、相位 $-\frac{\pi }{2}=-1.570796$（$\sin$ 写作 $\frac{e^{i\omega t}-e^{-i\omega t}}{2i}$，低频一侧的系数是 $-\frac{i}{2}$，相位 $-\frac{\pi }{2}$）；$\cos$ 项给出频率 $7$ 处幅度 $\frac{0.8}{2}=0.4$、相位 $0$；直流分量 $0.3$ 出现在 $k=0$。输出与这些值逐位一致，并且负频率一侧（$k=57$ 对应频率 $-7$，$k=61$ 对应频率 $-3$）的幅度与正频率相同、相位相反，正是 $X_{N-k}=\overline{X_k}$ 的表现。

实验四：正交展开的误差衰减

用 Chebyshev 部分和与 Legendre 最小二乘拟合逼近三个函数，比较误差衰减的形态。

import math
import numpy as np

def cheb_coeffs(fun, n):
    """用 Chebyshev-Gauss 节点计算 Chebyshev 展开系数"""
    x = np.cos(np.pi * (np.arange(2000) + 0.5) / 2000)
    w = np.pi / 2000
    a = np.zeros(n + 1)
    for i in range(n + 1):
        a[i] = (2.0 / math.pi) * float(np.sum(fun(x) * np.cos(i * np.arccos(x)) * w))
    a[0] /= 2.0
    return a

xfine = np.linspace(-0.999, 0.999, 2001)
funcs = [
    ("e^x", np.exp),
    ("|x|", np.abs),
    ("1/(1+25x^2)", lambda x: 1.0 / (1.0 + 25.0 * x * x)),
]
for name, fun in funcs:
    print(f"\nf(x) = {name}，Chebyshev 部分和的最大误差")
    print(f"{'n':>6} {'误差':>14}")
    for n in [2, 4, 8, 16, 32, 64]:
        a = cheb_coeffs(fun, n)
        err = float(np.max(np.abs(np.polynomial.chebyshev.chebval(xfine, a) - fun(xfine))))
        print(f"{n:>6} {err:>14.4e}")

print("Legendre 展开（均匀网格上的最小二乘拟合）")
for name, fun in funcs:
    print(f"  f(x) = {name}")
    for deg in [4, 8, 16, 32]:
        c = np.polynomial.legendre.legfit(xfine, fun(xfine), deg)
        err = float(np.max(np.abs(np.polynomial.legendre.legval(xfine, c) - fun(xfine))))
        print(f"    度 {deg:>3}: 最大误差 = {err:.4e}")

输出：

f(x) = e^x，Chebyshev 部分和的最大误差
     n             误差
     2     4.9901e-02
     4     5.8133e-04
     8     1.1575e-08
    16     1.7764e-15
    32     3.5527e-15
    64     1.5987e-14

f(x) = |x|，Chebyshev 部分和的最大误差
     n             误差
     2     2.1221e-01
     4     1.2732e-01
     8     7.0736e-02
    16     3.7449e-02
    32     1.9294e-02
    64     9.7984e-03

f(x) = 1/(1+25x^2)，Chebyshev 部分和的最大误差
     n             误差
     2     5.4027e-01
     4     3.6311e-01
     8     1.6401e-01
    16     3.3462e-02
    32     1.3929e-03
    64     2.4134e-06

Legendre 展开（均匀网格上的最小二乘拟合）
  f(x) = e^x
    度   4: 最大误差 = 1.1947e-03
    度   8: 最大误差 = 3.0510e-08
    度  16: 最大误差 = 5.8342e-14
    度  32: 最大误差 = 1.7392e-13
  f(x) = |x|
    度   4: 最大误差 = 1.1713e-01
    度   8: 最大误差 = 6.7257e-02
    度  16: 最大误差 = 3.6403e-02
    度  32: 最大误差 = 1.8998e-02
  f(x) = 1/(1+25x^2)
    度   4: 最大误差 = 3.3040e-01
    度   8: 最大误差 = 1.4944e-01
    度  16: 最大误差 = 3.0499e-02
    度  32: 最大误差 = 1.3637e-03

三组数据给出三种衰减形态，与光滑程度的阶梯一一对应。

$e^x$ 在复平面上解析，误差按几何速率下降：$4.99\times {10}^{-2}$、$5.81\times {10}^{-4}$、$1.16\times {10}^{-8}$，相邻比值的对数差近似为常数，$n$ 每加 $4$ 误差降约四个数量级。$n=16$ 达到 $1.78\times {10}^{-15}$（机器精度），之后不再改善（$n=32$ 的 $3.55\times {10}^{-15}$ 与 $n=64$ 的 $1.60\times {10}^{-14}$ 都在舍入噪声水平上，最后一行略高是因为部分和的计算累积了更多舍入误差）。这个函数用 $16$ 个基函数就能达到双精度极限，说明解析函数的谱展开效率极高。

$\left|x\right|$ 在 $x=0$ 处有拐角（一阶导数从 $-1$ 跳到 $1$），误差按代数速率下降：$2.12\times {10}^{-1}$、$1.27\times {10}^{-1}$、$7.07\times {10}^{-2}$、$3.74\times {10}^{-2}$、$1.93\times {10}^{-2}$、$9.80\times {10}^{-3}$。相邻比值是 $0.60$、$0.556$、$0.529$、$0.516$、$0.508$，逐渐稳定到 $0.5$，对应 $O(n^{-1})$。要降低一个数量级的误差需要把 $n$ 提高约十倍，与解析函数的情形形成对比。$\left|x\right|$ 的 Chebyshev 系数有精确表达（$\left|x\right|=\frac{2}{\pi }+\frac{4}{\pi }\sum _{j\ge 1}\frac{(-1{)}^j}{4j^2-1}{T}_{2j}(x)$，系数按 $\frac{1}{j^2}$ 衰减），截断误差是尾部求和，量级 $\frac{1}{n}$，与实测的比值 $0.5$ 相符。

$\frac{1}{1+25x^2}$ 解析但有距离 $0.2$ 的复奇点，误差按几何速率下降但速率取决于 $\rho =1.2198$：$5.40\times {10}^{-1}$、$3.63\times {10}^{-1}$、$1.64\times {10}^{-1}$、$3.35\times {10}^{-2}$、$1.39\times {10}^{-3}$、$2.41\times {10}^{-6}$。前几行的比值较大（$0.67$ 与 $0.45$），$n\ge 16$ 后渐近规律建立，比值降到 $0.042$ 与 $0.0017$（对应 ${\rho }^{-1}$ 的更高次幂在起作用：$n$ 每加 $16$ 是 ${\rho }^{-16}={1.2198}^{-16}=0.0417$，与 $n=16\to 32$ 的实测比值 $\frac{0.0013929}{0.033462}=0.0416$ 完全一致）。$n=64$ 时的理论量级 ${\rho }^{-64}=2.98\times {10}^{-6}$ 与实测的 $2.41\times {10}^{-6}$ 相符。

Legendre 一列给出同样的阶梯，只是常数因子不同：$e^x$ 在次数 $16$ 达到 $5.83\times {10}^{-14}$，$\left|x\right|$ 的误差按 $O(n^{-1})$ 下降（比值 $0.574$、$0.541$、$0.522$），$\frac{1}{1+25x^2}$ 从 $3.30\times {10}^{-1}$ 降到 $1.36\times {10}^{-3}$（次数 $4\to 32$，每加 $8$ 误差降到约 $\frac{1}{12}$，与 Chebyshev 展开的规律同源，因为收敛速率都由同一个复奇点决定）。两套基在最大误差上的差别通常在 $2$ 倍以内，选择时更多考虑别的因素：Chebyshev 基与插值节点、快速 Chebyshev 变换（类似 FFT 的算法）配套，Legendre 基与 Gauss-Legendre 求积、球谐函数展开配套。

5.7.7 本节小结

正交展开的三层内容

本节把幂级数与 Fourier 分析组织为三层。幂级数层给出 Cauchy-Hadamard 公式 $R=1/lim sup|a_k{|}^{1/k}$（实验一的三组系数分别给出半径 $\mathrm{\infty }$、$2$、$\frac{1}{e}=0.367879$，第三组的比值 ${(1+\frac{1}{k})}^k$ 以 $1-\frac{1}{2k}$ 的速率逼近 $e$，$k=100$ 时为 $2.704814$）、Abel 定理（内部闭子区间一致收敛）、逐项求导与逐项积分（半径不变，系数唯一决定和函数）以及解析性的判别（$e^{-1/x^2}$ 说明任意阶可导不等于解析，复奇点距离决定半径）。内积层建立函数空间的正交展开框架：$L^2$ 内积与均方距离、正交系的系数解耦（$c_k=\frac{⟨f,{\varphi }_k⟩}{⟨{\varphi }_k,{\varphi }_k⟩}$）、最佳均方逼近、Bessel 不等式（系数能量不超过函数能量）与完备性下的 Parseval 等式（能量守恒）。具体展开层给出两条实例。Fourier 级数把函数按三角系展开，三种收敛的含义不同：均方收敛对任何 $L^2$ 函数成立，逐点收敛在间断点处取左右极限的平均，一致收敛要求函数连续且周期延拓连续；Gibbs 现象的数据（方波最大过冲趋于 $1.178980=\frac{2}{\pi }\mathrm{Si}(\pi )$，相对跳跃 $8.949\mathrm{\%}$，与 $N$ 无关；同时均方误差按 $\sqrt{8/(\pi N)}$ 衰减，$N=501$ 时实测 $0.071333$）说明最大偏差与能量误差的收敛性质可以完全不同。DFT 与 FFT 把连续展开搬到有限数据上：DFT 是 $O(N^2)$ 的矩阵运算，FFT 用分治降到 $O(N\log N)$，实验三测得 $N=4096$ 时朴素 DFT 用 $0.406679$ 秒、FFT 用 $0.000029$ 秒，比值 $14023$，结果一致到 ${10}^{-10}$；频谱的共轭对称 $X_{N-k}=\overline{X_k}$ 与幅度公式在信号例子中逐位验证（$1.5\sin 3t+0.8\cos 7t+0.3$ 给出幅度 $0.75$、$0.4$ 与直流 $0.3$）。正交多项式层把插值节点（5.6 节）、Gauss 求积节点（5.5 节）与正交系的零点统一起来，并给出展开系数衰减与光滑程度的完整阶梯：解析函数几何衰减（$e^x$ 用 $16$ 项达机器精度）、含拐角函数按 $O(n^{-1})$（$\left|x\right|$ 的比值稳定在 $0.5$）、含复奇点函数按 ${\rho }^{-n}$ 几何衰减（$\frac{1}{1+25x^2}$ 的比值 $0.0416$ 与 ${\rho }^{-16}=0.0417$ 一致）。

正交展开处理的是「用振动或多项式分解函数」，下一节把微积分推广到多元。梯度、Hesse 矩阵、隐函数定理与 Lagrange 乘子构成了最优化与机器学习中大量方法的数学基础，其中矩阵求导的内容与第一章的矩阵微积分一节直接对应。

练习题

第 1 题 概念推导

(a) 求 $\sum _{k=0}^{\mathrm{\infty }}\frac{(x-1{)}^k}{3^k}$ 的收敛半径与收敛区间。(b) 用 Parseval 等式计算 $\sum _{k=1}^{\mathrm{\infty }}\frac{1}{k^2}$。(c) 说明为什么实解析比「任意阶可导」更强，并给出一个反例。

参考答案

(a) 令 $t=x-1$，级数为 $\sum \frac{t^k}{3^k}$，系数 $a_k=3^{-k}$，$|a_k{|}^{1/k}=\frac{1}{3}$，故 $R=3$。收敛区间为 $|x-1|<3$，即 $(-2,4)$。端点检验：$x=4$ 时级数为 $\sum 1^k=\sum 1$ 发散；$x=-2$ 时级数为 $\sum (-1{)}^k$，部分和在 $0$ 与 $1$ 之间跳动，发散。收敛区间是开区间 $(-2,4)$。在区间内和函数是 $\frac{1}{1-(x-1)/3}=\frac{3}{4-x}$。

(b) 取 $f(x)=x$ 在 $[-\pi ,\pi ]$ 上，其 Fourier 系数为 $b_k=\frac{2(-1{)}^{k+1}}{k}$，$a_k=0$。由 Parseval 等式

$$\parallel f{\parallel }_2^2={\int }_{-\pi }^{\pi }{x}^{2}dx=\frac{2{\pi }^3}{3}$$

而 $\pi \sum _{k=1}^{\mathrm{\infty }}{b}_k^2=\pi \sum \frac{4}{k^2}=4\pi \sum \frac{1}{k^2}$。两式相等得

$$\sum _{k=1}^{\mathrm{\infty }}\frac{1}{k^2}=\frac{2{\pi }^3/3}{4\pi }=\frac{{\pi }^2}{6}$$

这条推导是 Parseval 等式最经典的应用，也说明 Fourier 展开可以用来计算数项级数的和。取 $f(x)=x^2$ 并用同样方法可以得到 $\sum \frac{1}{k^4}=\frac{{\pi }^4}{90}$。

(c) 实解析要求 Taylor 级数在展开点的某个邻域内收敛到函数本身，这是一条关于余项的条件；任意阶可导只要求各阶导数存在。$f(x)=e^{-1/x^2}$（$x\ne 0$）与 $f(0)=0$ 满足 $f^{(k)}(0)=0$ 对一切 $k$ 成立（每一阶导数都是 $e^{-1/x^2}$ 乘以 $x$ 的有理函数，趋于零的速度快于任何幂函数），因此 Taylor 级数恒为零，与 $f$ 只在 $x=0$ 处相等。这个函数在 $x=0$ 处任意阶可导但不解析。对比复情形：复可导（全纯）蕴含解析，所以反例只能在实轴上构造。这条差别的实用含义是：用 Taylor 展开近似一个实函数时，需要单独确认收敛性，不能只检查「导数存在」。

第 2 题 计算推理

(a) 求方波 $f(x)=\mathrm{sgn}(x)$ 在 $[-\pi ,\pi ]$ 上的 Fourier 系数，并由此写出 $\frac{\pi }{4}$ 的级数表示。(b) 用 $\left|X_k\right|$ 与相位解释：为什么实数信号的 DFT 有一半是冗余的。(c) 估计用 Chebyshev 部分和逼近 $f(x)=e^x$ 在 $[-1,1]$ 上达到 ${10}^{-10}$ 精度所需的项数，并与 5.4 节的 Taylor 展开比较。

参考答案

(a) $f$ 是奇函数，$a_k=0$。正弦系数

$$b_k=\frac{1}{\pi }{\int }_{-\pi }^{\pi }\mathrm{sgn}(x)\sin kxdx=\frac{2}{\pi }{\int }_0^{\pi }\sin kxdx=\frac{2}{\pi }\cdot \frac{1-(-1{)}^k}{k}$$

$k$ 为偶数时 $b_k=0$，$k$ 为奇数时 $b_k=\frac{4}{\pi k}$，因此 $f(x)\sim \frac{4}{\pi }\sum _{k\text{奇}}\frac{\sin kx}{k}$。取 $x=\frac{\pi }{2}$：$f\left(\frac{\pi }{2}\right)=1$，而 $\sin \frac{k\pi }{2}$ 在 $k=1,5,9,\dots$ 时为 $1$，在 $k=3,7,11,\dots$ 时为 $-1$，因此

$$1=\frac{4}{\pi }(1-\frac{1}{3}+\frac{1}{5}-\frac{1}{7}+\cdots )⟹\frac{\pi }{4}=1-\frac{1}{3}+\frac{1}{5}-\frac{1}{7}+\cdots$$

这就是 Leibniz 级数。它是交错级数，收敛很慢（误差按 $\frac{1}{2N}$ 下降，要得到三位小数需要约 $500$ 项），与 Fourier 级数在间断点附近的慢收敛是同一个现象。

(b) 实数序列满足 $\overline{X_k}=\sum x_n{e}^{+2\pi ikn/N}=X_{N-k}$，因此幅度满足 $|X_{N-k}|=|X_k|$，相位满足 $\arg X_{N-k}=-\arg X_k$。长度为 $N$ 的实数序列有 $N$ 个自由度，DFT 的输出有 $N$ 个复数（$2N$ 个实数），两者之间的差正好是这条对称性给出的 $N$ 个约束，因此只有大约一半的频谱信息是独立的（$k=0$ 与 $k=\frac{N}{2}$ 例外，这两个频率的虚部为零）。实用含义：实信号的频谱存储与计算可以只用一半，快速算法中也有专门的实数 FFT（rfft）利用这条对称性把计算量再减半。

(c) $e^x$ 在 $[-1,1]$ 上的 Chebyshev 展开系数是 $c_k=2I_k(1)$（$I_k$ 是第一类修正 Bessel 函数，$k\ge 1$；$c_0=I_0(1)$），近似地 $c_k\approx \frac{1}{k!2^{k-1}}$（量级估计）。截断误差约等于被舍弃的首项，要求 $\frac{1}{n!2^{n-1}}\le {10}^{-10}$：$n=10$ 时 $\frac{1}{10!\cdot 2^9}=5.4\times {10}^{-9}$，$n=12$ 时 $\frac{1}{12!\cdot 2^{11}}=1.0\times {10}^{-11}$，因此约需 $12$ 项。实验四的数据显示 $n=16$ 时达到 $1.78\times {10}^{-15}$，与这个估计一致。作为对照，5.4 节的 Taylor 展开（Maclaurin 多项式）在 $x=1$ 处达到 ${10}^{-10}$ 需要 $n\approx 13$（误差 $\frac{e}{(n+1)!}$），两者的项数量级相同，但 Taylor 展开的误差在端点 $x=1$ 处最大而在 $x=0$ 附近极小（误差分布很不均匀），Chebyshev 展开的误差在整段上接近等幅。整段精度要求高时 Chebyshev 展开更合适，只需要端点附近的局部精度时 Taylor 展开更简单。

第 3 题 代码验证

用本节代码做三组实验并解释：(a) 在实验二里把方波换成周期为 $1$、在 $[0,\frac{1}{2})$ 上取 $1$、在 $[\frac{1}{2},1)$ 上取 $-1$ 的方波，检查过冲的极限值是否改变；(b) 在实验三里把输入信号的频率换成非整数（例如 $3.4$）并观察频谱；(c) 在实验四里把 $\left|x\right|$ 换成 ${\left|x\right|}^3$（二阶连续可导），观察 Chebyshev 误差的衰减速率变化。

参考答案

(a) 过冲的极限值不变。Gibbs 现象是局部的：跳跃处附近的振荡形状由截断的高频分量决定，与跳跃的位置、幅度以及函数在别处的取值无关。改变周期只改变 Fourier 系数的具体形式（频率间隔变成整数倍的 $\frac{1}{T}$），跳跃幅度不变则过冲的绝对幅度不变（$0.178980\times$ 跳跃幅度的一半），相对比例仍是 $8.949\mathrm{\%}$。数值上表现为：最大过冲值仍是 $1.178980$，只是过冲出现的位置随跳跃点位置平移。这条结论的实用含义是，只要函数存在跳跃，用有限的 Fourier 项逼近就一定出现过冲，改变展开区间或周期不能解决。

(b) 频率 $3.4$ 不在 DFT 的频率格点上（格点是 $\frac{2\pi k}{N}$ 对应的整数 $k$），信号的能量会泄漏到邻近的频率分量上：$k=3$ 与 $k=4$ 附近都出现较大幅度，且没有一条谱线等于真实幅度 $0.75$，同时幅度随 $k$ 的变化呈现由 $\mathrm{Dirichlet}$ 核（即 $\frac{\sin (N\theta /2)}{\sin (\theta /2)}$）决定的旁瓣结构。这是频谱泄漏现象：DFT 隐含地把信号视为周期延拓，非整周期信号的周期延拓产生跳跃，跳跃的能量按 $O(1/|k-k_0|)$ 的速度扩散到整个频谱。处理方法是加窗（用汉宁窗、汉明窗等平滑因子压低旁瓣）、增加采样长度（泄漏的旁瓣幅度随 $N$ 减小）或用零填充（改善谱线的插值密度，但不减少泄漏）。这条现象也是短时傅里叶分析要仔细选择窗长的原因。

(c) ${\left|x\right|}^3$ 在 $x=0$ 处有二阶导数连续（$f^{″}=6\left|x\right|$，连续但不可导），按系数的衰减阶梯应当落在 $\left|x\right|$ 与解析函数之间。Chebyshev 系数由函数的光滑程度决定：${\left|x\right|}^3$ 的二阶导数连续，对应系数按 $O(k^{-4})$ 衰减（$p$ 阶连续可导给出 $O(k^{-p-1})$ 量级），截断误差按 $O(n^{-3})$ 下降。数值上的检验方法是观察「$n$ 翻一倍时误差降到多少分之一」：$\left|x\right|$ 的比值趋于 $\frac{1}{2}$（$O(n^{-1})$），${\left|x\right|}^3$ 的比值应当趋于 $\frac{1}{8}$（$O(n^{-3})$）。这条练习给出一种从数据反推光滑程度的方法，实际应用中用来判断「继续增加基函数是否值得」。

第 4 题 综合应用

某音频处理程序需要对长度 ${10}^6$ 的采样序列做以下三步运算：低通滤波、计算与另一段信号的互相关、计算功率谱。请比较三种实现方案并给出推荐：直接按定义计算、用 FFT 加频域乘法、用分段处理（overlap-add）。同时说明各自的复杂度与适用条件。

参考答案

方案一：按定义计算。低通滤波是卷积运算，直接实现需要 $O(N\cdot M)$ 次乘加（$M$ 为滤波器长度）；互相关按定义需要 $O(N^2)$ 次运算；功率谱需要 $O(N^2)$ 次（朴素 DFT）。$N={10}^6$ 时互相关与功率谱的运算量在 ${10}^{12}$ 量级，秒级到分钟级，不可接受。方案一只在小 $N$ 或滤波器极短（$M$ 不超过几十）时可用，此时直接卷积的常数因子比 FFT 小，反而更快。

方案二：FFT 加频域乘法。三条运算都可以用卷积定理统一处理：滤波与互相关是卷积的形式，通过 $O(N\log N)$ 的 FFT 转到频域后做逐点乘法即可；功率谱本身就是 $|X_k{|}^2$，一次 FFT 足够。$N={10}^6$ 时 $N\log N\approx 2\times {10}^7$ 次运算，在毫秒级完成（本节实验三的数据支持这个量级：$N=4096$ 时 FFT 用 $0.000029$ 秒，按 $N\log N$ 外推 $N={10}^6$ 时约 $0.02$ 秒）。需要注意的细节是循环卷积与线性卷积的差别：DFT 给出的是循环卷积，要得到线性卷积需要把序列零填充到长度 $N+M-1$，否则尾部会绕回头部造成混叠。

方案三：分段处理（overlap-add 或 overlap-save）。把长序列切成若干段（每段长度 $L$，与滤波器长度 $M$ 匹配），每段用方案二处理，再把结果按重叠部分相加。优点是内存占用低（不需要同时存下整个 ${10}^6$ 长度的频谱，每段只需 $O(L)$ 内存）与延迟低（可以流式处理，边采边算）；代价是段边界处的处理需要额外的重叠计算，总运算量比一次性 FFT 略高（约高 $\frac{L}{L-M}$ 倍）。适用于实时音频处理、滤波器很长或内存受限的场合。

推荐：${10}^6$ 长度、可离线处理时用方案二，实现最简单，一次 FFT 解决三个需求（滤波与互相关的频域乘法、功率谱的模平方），总代价 $O(N\log N)$。要求实时或内存受限时用方案三，分段长度取 $L=4096$ 到 $16384$（$2$ 的幂便于 FFT，同时与常见音频缓冲区大小兼容）。三条运算共用的实现要点是复用 FFT 结果：功率谱与互相关都来自同一个变换，预先算一次频谱即可，避免重复变换。报告中应当写出所用的窗函数（若有）、零填充长度与分段策略，因为这些选择直接影响频谱泄漏与边界效应。

常见错误

错误 1 · 用解析公式的和函数判断幂级数在半径外的行为

原因：解析公式 $\frac{1}{1-x}$ 在 $x\ne 1$ 处都有定义，而级数 $\sum x^k$ 只在 $|x|<1$ 内收敛。实验一的数据把差别写得直接：$x=1.01$ 时部分和（$N=100$）是 $173.18619677$，而解析公式给出 $-100$；$x=1.1$ 时部分和是 $151576.73573805$，解析公式给出 $-10$。两者符号相反、量级无关。把解析延拓当作级数的和，会在半径外得到完全错误的结果。另一类相关的错误是把端点 $|x|=R$ 当作「一定收敛」或「一定发散」：端点需要单独判断（$\sum \frac{x^k}{k}$ 在 $x=-1$ 收敛、在 $x=1$ 发散）。

解决：使用幂级数前先算收敛半径（Cauchy-Hadamard 或比值法），并检验所用 $x$ 是否在半径内；需要半径外的信息时用解析延拓的公式，但要清楚它与级数求和是两件事。端点情形单独代入检验（比较判别法、Leibniz 判别法或直接看通项）。

错误 2 · 认为增大 Fourier 部分和的项数可以消除过冲

原因：Gibbs 现象的过冲幅度趋于非零极限 $\frac{2}{\pi }\mathrm{Si}(\pi )-1=0.178980$，与项数无关。实验二的数据从 $N=5$ 到 $N=501$，过冲从 $1.188357$ 只降到 $1.178980$，$100$ 倍的数据量几乎没有改变这个值（变化量 $0.009$）。过冲的位置随 $N$ 移向间断点（$x\approx \frac{\pi }{N}$），宽度趋于零，但最大偏差不趋于零。若在信号重建后看到固定的振铃，增加项数不会消除它。

解决：认识到过冲的根源是间断（Fourier 级数在跳跃处不可能一致收敛），处理方式是改变展开或改变函数：加窗压低高频系数的幅度（用 ${\sigma }_k=\frac{\sin (k\pi /N)}{k\pi /N}$ 一类的因子乘在系数上，代价是分辨率下降）、对信号做平滑或总变差型正则化、或者改用不产生过冲的基（例如用 $L^1$ 型基或分段多项式）。报告中说明所用的处理方法与它带来的分辨率损失。

错误 3 · 把循环卷积当作线性卷积

原因：DFT 的卷积定理给出的是循环卷积：$\mathrm{IDFT}(X_k{Y}_k)$ 等于把两个序列视为周期延拓后的卷积，长度为 $min(N,M)$ 的尾端会绕回到开头叠加。要求线性卷积（普通卷积）时，若直接把两个长度 $N$、$M$ 的序列做 $N$ 点 FFT 相乘再逆变换，得到的结果在尾部有混叠误差（能量从结尾绕到开头）。这个错误在滤波器实现中表现为开头出现「回声」或数值异常，很容易被误判为浮点问题。实验三只验证了 FFT 与朴素 DFT 的数值一致性（偏差 ${10}^{-12}$）与频谱的对称性，卷积定理的使用需要额外的零填充步骤。

解决：计算线性卷积时把两个序列零填充到长度 $N+M-1$（或下一个 $2$ 的幂以便 FFT），再做逐点乘法。要处理很长的序列时用 overlap-add 或 overlap-save 分段法，每段内部仍按 $L+M-1$ 的长度做变换。实现后用「直接卷积与 FFT 卷积的结果逐位对比」做验证，这是最省事的检查方式。

错误 4 · 把展开系数的衰减速率当作与函数无关的性质

原因：系数衰减速率由函数的光滑程度决定，同一套基对不同函数的效率差别很大。实验四的 Chebyshev 展开：$e^x$（解析）用 $16$ 项达到 ${10}^{-15}$，$\left|x\right|$（有拐角）用 $64$ 项只达到 ${10}^{-2}$，$\frac{1}{1+25x^2}$（有距离 $0.2$ 的复奇点）用 $64$ 项达到 ${10}^{-6}$。三者的项数需求相差十几个数量级。若按「解析函数用 16 项」的经验去处理含拐角的函数，会得到远低于预期的精度，或者误判为程序有错。

解决：先判断函数的光滑程度（有几阶连续导数、有没有间断、复奇点离区间多远），再据此估计收敛速率与所需项数。实用的检验方法是计算连续几档截断误差的比值：比值稳定在 $\frac{1}{2}$ 说明 $O(n^{-1})$，稳定在 $\frac{1}{4}$ 说明 $O(n^{-2})$，稳定在某个小于 $1$ 的常数说明几何收敛。若比值趋于 $1$，说明当前的基与函数不匹配（例如对有尖峰的周期信号用少量 Fourier 项），应当更换基或先做变量替换去掉奇性。