九点差分格式截断误差如何推导？

Question

在使用九点差分格式求解二维泊松方程时，如何系统推导其截断误差？常见问题在于：许多教材仅给出五点格式的误差分析，而九点格式涉及更多邻点组合，需利用泰勒展开将中心点周围八个节点的函数值展开至高阶项，进而代入差分公式。关键难点是如何合并同类项并识别主导误差项，尤其当格式通过组合拉普拉斯算子与双调和项提高精度时，如何确认其截断误差达到 $O(h^4)$？这一推导过程对理解格式精度提升机制至关重要。

Answer 1

使用九点差分格式求解二维泊松方程的截断误差系统推导

1. 问题背景与动机：从五点到九点格式的精度跃迁

在数值求解偏微分方程（PDE）中，二维泊松方程：

\[ \nabla^2 u = f(x, y) \]

是常见的模型问题。传统五点差分格式具有二阶精度 $O(h^2)$，其截断误差分析在多数教材中已有详述。然而，当需要更高精度时，九点差分格式因其潜在的四阶收敛性而备受关注。

九点格式通过引入更多邻点信息（包括对角点），构造出包含拉普拉斯算子和双调和项的组合形式，从而提升局部逼近阶数。但其截断误差的系统推导较为复杂，涉及多变量泰勒展开、高阶导数组合与对称性利用。

2. 九点差分格式的基本形式

设网格步长为 $h$，中心点为 $(x_0, y_0)$，其八个邻点构成如下模板：

$u_{-1,-1}$	$u_{0,-1}$	$u_{1,-1}$
$u_{-1,0}$	$u_{0,0}$	$u_{1,0}$
$u_{-1,1}$	$u_{0,1}$	$u_{1,1}$

九点格式的一般形式可写为：

\[ \frac{1}{6h^2} \left( 4(u_{1,0} + u_{-1,0} + u_{0,1} + u_{0,-1}) + (u_{1,1} + u_{-1,1} + u_{1,-1} + u_{-1,-1}) - 20u_{0,0} \right) = f_{0,0} + T.E. \]

其中 $T.E.$ 表示截断误差，目标是确定其阶数。

3. 泰勒展开：构建局部逼近的基础工具

将每个邻点函数值在中心点 $(x_0, y_0)$ 处进行二维泰勒展开至六阶项：

• $u_{\pm1,0} = u \pm h u_x + \frac{h^2}{2} u_{xx} \pm \frac{h^3}{6} u_{xxx} + \frac{h^4}{24} u_{xxxx} \pm \cdots$
• $u_{0,\pm1} = u \pm h u_y + \frac{h^2}{2} u_{yy} \pm \frac{h^3}{6} u_{yyy} + \frac{h^4}{24} u_{yyyy} \pm \cdots$
• $u_{\pm1,\pm1} = u \pm h(u_x + u_y) + \frac{h^2}{2}(u_{xx} + 2u_{xy} + u_{yy}) \pm \frac{h^3}{6}(u_{xxx} + 3u_{xxy} + 3u_{xyy} + u_{yyy}) + \cdots$

共需展开9个点，保留至 $h^6$ 阶以确保能识别 $O(h^4)$ 截断误差。

4. 代入差分算子并合并同类项

将所有展开式代入九点格式左侧，逐项相加：

令 S = 4(u_{1,0} + u_{-1,0} + u_{0,1} + u_{0,-1}) 
       + (u_{1,1} + u_{-1,1} + u_{1,-1} + u_{-1,-1}) 
       - 20u

计算各项系数：

导数项	来自轴向点贡献	来自对角点贡献	总系数
$u_{xx}, u_{yy}$	$4 \cdot \frac{h^2}{2} \cdot 2 = 4h^2$	$\frac{h^2}{2} \cdot 4 = 2h^2$	$6h^2(u_{xx}+u_{yy})$
$u_{xxxx}, u_{yyyy}$	$4 \cdot \frac{h^4}{24} \cdot 2 = \frac{h^4}{3}$	$\frac{h^4}{24} \cdot 4 = \frac{h^4}{6}$	$\frac{h^4}{2}(u_{xxxx}+u_{yyyy})$
$u_{xxyy}$	0	$\frac{h^4}{2} \cdot 4 \cdot 2u_{xxyy}/2 = 2h^4 u_{xxyy}$	$2h^4 u_{xxyy}$

注意：由于对称性，所有奇数阶导数项（如 $u_x, u_{xxx}$）在求和后完全抵消。

5. 构造离散拉普拉斯与双调和项的关系

观察发现，九点格式本质上是以下组合：

\[ \Delta_h^{(9)} u = \frac{1}{6h^2} S = \nabla^2 u + \frac{h^2}{12} \nabla^4 u + O(h^4) \]

其中 $\nabla^4 u = \Delta(\Delta u) = u_{xxxx} + 2u_{xxyy} + u_{yyyy}$ 为双调和算子。

因此，原始泊松方程 $\nabla^2 u = f$ 的离散化变为：

\[ \Delta_h^{(9)} u = f + \frac{h^2}{12} \nabla^2 f + O(h^4) \]

这表明若 $f$ 光滑，可通过修正源项消除 $h^2$ 项，实现整体四阶精度。

6. 截断误差主导项识别与阶数确认

由上述展开可得截断误差表达式：

\[ T.E. = \Delta_h^{(9)} u - \nabla^2 u = \frac{h^2}{12} \nabla^4 u + O(h^4) \]

但进一步分析高阶项发现，由于格式设计巧妙利用了对称性，$h^2$ 项被精确补偿，实际误差始于 $h^4$ 阶。

关键在于：九点格式并非简单逼近 $\nabla^2 u$，而是构造了一个满足：

\[ \Delta_h^{(9)} u = \nabla^2 u + C h^4 \mathcal{L}_6 u + O(h^6) \]

的算子，其中 $\mathcal{L}_6$ 是六阶微分算子。因此，在适当条件下，其局部截断误差为 $O(h^4)$。

7. 常见技术问题与解决方案

1. 泰勒展开遗漏交叉导数项：应完整展开至四阶以上，特别注意 $u_{xxyy}, u_{xxxy}$ 等混合项。
2. 误判主导误差阶数：仅看前几项可能误认为有 $h^2$ 项，需验证其系数是否因格式系数选择而恰好为零。
3. 边界处理破坏精度：即使内部达到 $O(h^4)$，边界若用低阶格式会拖累整体收敛率。
4. 非均匀网格适应性差：九点格式在非结构网格上难以直接推广，需重构权重。

解决方案包括：采用符号计算工具（如SymPy）辅助展开；使用Richardson外推验证实际收敛阶；结合紧致差分思想设计边界格式。

8. 可视化推导流程：mermaid 流程图展示

graph TD
    A[开始] --> B[定义九点模板]
    B --> C[对9个点做泰勒展开至h^6]
    C --> D[代入差分公式S]
    D --> E[按导数阶数合并同类项]
    E --> F[识别u_xx, u_yy项→形成∇²u]
    F --> G[提取h²项: 检查是否抵消]
    G --> H[分析h⁴项: ∇⁴u或更高]
    H --> I[得出截断误差表达式]
    I --> J[结论: T.E. = O(h⁴)]

9. 实际应用中的精度验证方法

理论推导需配合数值实验验证。常用方法：

• 选取解析解 $u(x,y)$，计算对应 $f = \nabla^2 u$
• 在不同网格 $h, h/2, h/4$ 上求解
• 计算最大误差 $E_h = \max |u_h - u_{exact}|$
• 拟合 $\log E_h \sim p \log h$，斜率 $p$ 应接近4

例如取 $u = \sin(\pi x)\sin(\pi y)$，则 $f = -2\pi^2 u$，易于实现。

10. 扩展思考：高阶格式的设计哲学

九点格式的成功启示我们：

1. 增加模板点数可引入更高阶微分信息
2. 通过对称性消除低阶误差项
3. 组合基本算子（如$\nabla^2$与$\nabla^4$）可构造高阶逼近
4. 精度提升依赖于解的光滑性假设

这一思想延伸至紧致差分、谱方法、有限体积高阶重构等领域。