12、热方程的数值解法

热方程的数值解法

1. 隐式方法

隐式方法基于向后时间差分方程 $\frac{u^{\ell + 1}-u^{\ell}}{dt}=cD^{2}u^{\ell + 1}$，每一个时间步都需要求解一个系统 $B_{d}u^{\ell + 1}=u^{\ell}$，其中 $B_{d}=I - dD^{2}$。以下是实现该方法一个时间步的代码：

u = B \ u;

这个方法是稳定的，较小的 $\alpha$ 值可能会得到更好的结果，但较大的值也不会导致灾难性的失败。不过，历史上求解系统的实现较为困难且执行速度慢，现在使用现代计算机上的求解器（如 \ ），该方法也变得同样“容易”。

示例 4.4 ：设 $d = \frac{2}{2+\sqrt{2}}$，考虑 $(0, 1)$ 上的 $n = 4$ 点网格，边界条件为 0 - D。初始温度 $u_{0}=(\frac{\sqrt{2}}{2}, 1, \frac{\sqrt{2}}{2})$，求解隐式方法的一个时间步，并将结果振幅与 $t = \Delta t$ 时的精确已知振幅以及 $B_{d}^{-1}$ 的特征值 $\lambda_{B_{d}^{-1},1}$ 进行比较。

Bd = (1/(2 + sqrt(2)))*[6 + sqrt(2) -2 -2; -2 6 + sqrt(2) -2; -2 -2 6 + sqrt(2)];
u0 = [sqrt(2)/2; 1; sqrt