2. 有限差分基础与模型方程

稳定性分析不是凭空展开的，它必须建立在一套明确的离散框架之上。本节的任务是把“要分析什么”这件事说清楚：从导数离散、网格划分、截断误差与相容性这些基础概念出发，以热传导方程和波动方程为例，建立半离散与全离散格式，并明确稳定性分析所针对的对象——即离散格式在时间推进过程中误差是否被放大。

一、从导数到差分：网格与离散算子

有限差分的第一步是把连续区域替换为离散网格。以一维区间 [0, L] 为例，取空间步长 Δx = L/N，节点 x_j = jΔx（j = 0,1,…,N）；时间方向取步长 Δt，时刻 tⁿ = nΔt。记 u_jⁿ ≈ u(x_j, tⁿ)。所有差分格式都在这张网格上定义，网格比 λ = Δt/Δx 或 r = Δt/Δx² 往往直接决定稳定性结论。

导数离散来自 Taylor 展开。对 u(x_j±Δx) 在 x_j 处展开并组合，可得三种常用空间算子：

• 一阶中心差分：(u_j+1 − u_j−1)/(2Δx) = u′ + O(Δx²)，二阶精度；
• 二阶中心差分：(u_j+1 − 2u_j + u_j−1)/Δx² = u″ + O(Δx²)，这是热传导方程空间离散的核心；
• 一阶向前/向后差分：(u_j+1 − u_j)/Δx 或 (u_j − u_j−1)/Δx，一阶精度，常用于时间方向或迎风处理。

时间方向同理：向前差分 u_t ≈ (u_jⁿ⁺¹ − u_jⁿ)/Δt 为一阶精度，中心差分 (u_jⁿ⁺¹ − u_jⁿ⁻¹)/(2Δt) 为二阶精度。选择哪一组算子，就决定了后面要分析的格式形态。

二、截断误差与相容性

把精确解代入差分格式，得到的余项称为截断误差。例如对热传导方程 u_t = αu_xx 的显式格式，截断误差为 O(Δt) + O(Δx²)。若 Δt、Δx → 0 时截断误差趋于零，则称格式与微分方程相容。相容性回答的是“格式逼近的是不是原方程”，它只涉及局部近似精度，与误差是否随时间增长无关。这一点很关键：一个相容的格式完全可能是不稳定的，因此相容性不能替代稳定性分析。

三、半离散与全离散：以两个模型方程为例

热传导方程 u_t = αu_xx。只离散空间、保留时间连续，得到半离散方程组：

du_j/dt = α (u_{j+1} − 2u_j + u_{j−1}) / Δx²,  j = 1,…,N−1

这是一组常微分方程，可写成 du/dt = A u，矩阵 A 为三对角、对称负定（配合边界条件），其特征值决定了半离散系统的行为。进一步对时间用向前 Euler，得到全离散显式格式：

u_j^{n+1} = u_j^n + r (u_{j+1}^n − 2u_j^n + u_{j−1}^n),  r = αΔt/Δx²

若时间改用向后 Euler，则得到全隐式格式 u_jⁿ⁺¹ − r(u_j+1ⁿ⁺¹ − 2u_jⁿ⁺¹ + u_j−1ⁿ⁺¹) = u_jⁿ，每步需解三对角线性方程组。Crank–Nicolson 则对时间取平均，二阶精度、无条件稳定，但需注意其在高频模态上的振荡特性。

波动方程 u_tt = c²u_xx。空间二阶中心差分后，时间用二阶中心差分，得到经典显式格式：

u_j^{n+1} = 2u_j^n − u_j^{n−1} + C² (u_{j+1}^n − 2u_j^n + u_{j−1}^n),  C = cΔt/Δx

这里 C 称为 Courant 数。该格式的稳定性条件为 C ≤ 1，即 Courant–Friedrichs–Lewy 条件，物理含义是数值依赖域必须覆盖解析依赖域。与热传导方程不同，波动方程是双曲型的，其稳定性约束来自“信息传播速度”而非“扩散耗散”。

四、稳定性分析的对象

综合上述例子，稳定性分析的对象可以精确表述为：离散格式的误差传播算子。设 Uⁿ 为数值解向量，uⁿ 为精确解在网格上的限制，误差 eⁿ = Uⁿ − uⁿ。对线性常系数格式，误差满足与数值解相同的齐次递推 eⁿ⁺¹ = G eⁿ（单步格式）或相应的多步关系，G 即放大矩阵。稳定性要求存在与 n 无关的常数 K，使 ‖Gⁿ‖ ≤ K 对所有 n 成立。

对半离散系统 du/dt = Au，稳定性对应 A 的特征值实部不超过零（或落在稳定域内）；对全离散格式，则要求 Δt 与 Δx 满足特定约束，使 G 的谱半径不超过 1（或一致有界）。这正是后续 von Neumann 分析、矩阵方法和能量方法要处理的核心问题。

五、一个可操作的判断流程

面对一个具体问题，可按以下顺序落实本节内容。以“α = 0.01、Δx = 0.05、要求显式格式稳定”为例：先由热传导显式格式条件 r = αΔt/Δx² ≤ 1/2 反解 Δt ≤ 0.5Δx²/α = 0.5×0.0025/0.01 = 0.125。若取 Δt = 0.1，则 r = 0.4，满足条件；若取 Δt = 0.2，则 r = 0.8 > 0.5，格式将出现非物理振荡并迅速发散。对波动方程，若 c = 1、Δx = 0.05，则要求 Δt ≤ Δx/c = 0.05，即 C ≤ 1。这些数值均为教学示例，用于说明参数与稳定性条件的对应关系。

需要强调的是，上述条件只对线性常系数、规则网格、特定格式成立。变系数、非均匀网格、边界处理和高维情形都需要重新分析，边界格式甚至可能比内点格式更严格。因此，明确“分析对象是哪一个放大矩阵、在什么范数下、对哪些模态成立”，是进入后续 von Neumann 分析与矩阵稳定性讨论的前提。