一维有限元残差后验估计

解出有限元矩阵方程之后，还有一个问题没有回答：这张网格上的折线，距离原来的连续解究竟多远？矩阵残差为零，只说明有限个测试方向上的平衡已经满足；它没有检查单元内部的弯曲，也没有消除折线在节点处的斜率变化。

残差后验估计从已经得到的离散解和已知荷载出发，为这些缺陷赋予正确的尺度，再把它们变成误差界。与需要未知解正则性范数的先验估计相比，它更适合回答一次具体计算是否足够精细。本条在一个可以把每一步证明到底的一维模型中建立这座桥；自适应有限元公理库自适应有限元的估计、标记与加密Adaptive finite element method · AFEM · Dörfler marking在局部荷载的一维 Poisson 问题上完成求解、残差估计、平方指标标记与局部二分，并比较相同误差和相同预算下的网格选择。随后用局部指标决定在哪里增加自由度。

形式陈述 ​

固定方程、空间和误差 ​

考虑

$$-u^{″}=f\text{于 }(0,1),u(0)=u(1)=0,f\in L^2(0,1).$$

取任意严格递增划分 $0=x_0<\cdots <x_N=1$，单元 $K_i=(x_{i-1},x_i)$ 的长度为 $h_i$。$V_h$ 是连续分片一次、端点为零的空间。按有限元方法公理库有限元方法Finite element method · FEM · Finite element assembly从弱形式与帽函数出发，完整构造局部单元、组装并求解小系统，再用能量与逼近解释误差阶和边界处理。精确计算积分，并精确求解离散系统，所得 $u_h\in V_h$ 满足

$${\int }_0^1{u}_h^{\prime }{v}_h^{\prime }dx={\int }_0^{1}f{v}_{h}dx(v_h\in V_h).$$

这里讨论的误差是能量误差

$$e=u-u_h,E=\Vert e^{\prime }{\Vert }_{L^2(0,1)},E_K=\Vert e^{\prime }{\Vert }_{L^2(K)}.$$

零端点条件使导数半范数成为 $H_0^1$ 上的范数。数值积分误差和未完成线性求解造成的误差暂不包括在 $E$ 中；后文会说明怎样区分代数误差。

单元内部与节点处的残差 ​

对任意 $v\in H_0^1(0,1)$，定义连续残差泛函

$$R(v)={\int }_0^{1}fvdx-{\int }_0^1{u}_h^{\prime }{v}^{\prime }dx.$$

设 $U_i=u_h(x_i)$，单元斜率和内部节点跳跃分别为

$$s_i=\frac{U_i-U_{i-1}}{h_i},J_i=s_{i+1}-s_i(1\le i<N).$$

逐单元积分得到

$$R(v)=\sum _K{\int }_{K}fvdx+\sum _{i=1}^{N-1}{J}_{i}v(x_i).$$

符号来自相邻端点贡献：左单元给出 $-s_{i}v(x_i)$，右单元给出 $+s_{i+1}v(x_i)$。单元内 $u_h^{″}=0$，所以强残差 $f+u_h^{″}$ 就是 $f$；但分布意义下的二阶导数还包含节点贡献，不能因为单元内二阶导数为零就略去跳跃。两端的测试函数为零，因此此处没有 Dirichlet 端点跳跃项。

取 $v=v_h$ 时 $R(v_h)=0$，并不要求上式中的每一项分别为零。对内部帽函数 ${\varphi }_i$，更精确的平衡是

$$J_i=-{\int }_0^{1}f{\varphi }_{i}dx.$$

这也说明 $J_i$ 是数值解斜率的跳跃，并非荷载 $f$ 本身的跳跃。

本条采用的权重与局部分配 ​

定义

$$V^2=\sum _K{h}_K^2\Vert f{\Vert }_{L^2(K)}^2,w_i=min(h_i,h_{i+1}),{\eta }^2=V^2+\sum _{i=1}^{N-1}{w}_i{J}_i^2.$$

将一个内部节点的跳跃贡献平均分给左右单元：

$${\eta }_{K_i}^2=h_i^2\Vert f{\Vert }_{L^2(K_i)}^2+\frac{1}{2}{w}_{i-1}{J}_{i-1}^2+\frac{1}{2}{w}_i{J}_i^2,$$

其中不存在的端点项直接省略。于是 ${\eta }^2=\sum _K{\eta }_K^2$，每个 ${\eta }_K$ 都可由当前解、荷载和网格计算。

这是本条明确选择的一维指标，权重使用相邻长度的较小值，便于在高度不均匀的网格上给出统一常数。多维残差估计也组合单元残差和界面通量跳跃，但界面测度、权重与局部邻域需要重新定义；这里的公式不是把高维面面积机械代成零而得到的特例。[1]

为什么它能控制误差 ​

一个一维特性：离散解等于节点插值 ​

$u\in H^2(0,1)$，因此可以定义节点插值 $I_{h}u$。任意 $v_h\in V_h$ 的导数在各单元上为常数，而

$${\int }_K(u-I_{h}u{)}^{\prime }dx=0.$$

所以 $a(u-I_{h}u,v_h)=0$。这与Galerkin 正交公理库Galerkin 方法Galerkin method · Galerkin discretization在有限维试探与检验空间中离散连续变分问题，并用 Galerkin 正交、Céa 准最优性及稳定条件组织误差。完全相同；离散解唯一，故 $u_h=I_{h}u$。它不是说整个函数已经精确，而是说 $e$ 在每个单元的两个端点都为零。在各单元内部，仍有 $-e^{″}=f$。

这个结论依赖一维、单位扩散系数和协调一次元。它是下述简洁证明的关键，不能直接移植到一般变系数或高维网格。

可靠性：估计量不会漏掉大误差 ​

长度为 $h_K$ 的区间上，零端点函数满足 Poincaré 不等式 $\Vert e{\Vert }_K\le (h_K/\pi )\Vert e^{\prime }{\Vert }_K$。局部分部积分后，

$$E_K^2={\int }_{K}fedx\le \Vert f{\Vert }_K\frac{h_K}{\pi }{E}_K.$$

约去非零的 $E_K$，再平方求和，得到

$$\boxed{{\displaystyle E\le \frac{V}{\pi }\le \frac{\eta }{\pi }.}}$$

这就是可靠性公理库残差、误差估计与停止准则Residual and error estimation · Stopping criterion区分可计算残差与未知真误差，并说明把缺陷转成误差界和停止证书所需的条件。：只要计算假设成立，$\eta /\pi$ 是能量误差的上界。对这个特殊模型，仅体积项 $V$ 就足以给出上界；保留跳跃项有助于呈现连续残差的完整结构，并为局部标记提供另一种分配信息。

荷载在每个单元上常值时：可以精确比较 ​

若 $f{|}_K=c_K$，写局部坐标 $t=x-x_{i-1}$。由 $-e^{″}=c_K$ 与零端点条件，

$$e(t)=\frac{c_K}{2}t(h_K-t),E_K^2={\int }_0^{h_K}{c}_K^2(h_K/2-t{)}^{2}dt=\frac{c_K^2{h}_K^3}{12}.$$

因此 $V^2=12E^2$。帽函数平衡进一步给出

$$J_i=-\frac{c_i{h}_i+c_{i+1}{h}_{i+1}}{2}.$$

由 $(a+b{)}^2\le 2(a^2+b^2)$ 以及 $w_i\le h_i,h_{i+1}$，每个节点贡献满足

$$w_i{J}_i^2\le \frac{1}{2}(c_i^2{h}_i^3+c_{i+1}^2{h}_{i+1}^3).$$

一个单元至多被相邻两个节点计入，故 $\sum _i{w}_i{J}_i^2\le V^2$。于是

$$\boxed{{\displaystyle \sqrt{12}E\le \eta \le \sqrt{24}E.}}$$

这里“荷载分片常值”必须相对于当前网格成立。若荷载断点落在单元内部，不能直接套用此精确公式。满足条件时，二分一个单元会把该单元的误差平方从 $c_K^2{h}_K^3/12$ 变成原来的 $1/4$；其他单元的误差不变。这为一次局部加密提供了可以独立核验的真误差答案。

数据振荡为何不能忽略 ​

把荷载分成平均值与未分辨部分 ​

对一般 $f\in L^2$，定义

$${\overline{f}}_K=\frac{1}{h_K}{\int }_{K}fdx,{\mathrm{osc}}^2=\sum _K{h}_K^2\Vert f-{\overline{f}}_K{\Vert }_K^2,{\overline{V}}^2=\sum _K{h}_K^2\Vert {\overline{f}}_K{\Vert }_K^2.$$

正交分解给出 $V^2={\overline{V}}^2+{\mathrm{osc}}^2$。单元上再把误差分成 $e=e_0+e_d$，其中两个函数都在单元端点为零，分别满足 $-e_0^{″}={\overline{f}}_K$ 和 $-e_d^{″}=f-{\overline{f}}_K$。前一项能精确积分，后一项应用刚才的可靠性证明，因此

$$\Vert e_0^{\prime }\Vert =\frac{\overline{V}}{\sqrt{12}},\Vert e_d^{\prime }\Vert \le \frac{\mathrm{osc}}{\pi }.$$

这里范数对全部单元求和。由三角不等式及其反向形式，

$$max(0,\frac{\overline{V}}{\sqrt{12}}-\frac{\mathrm{osc}}{\pi })\le E\le \frac{\overline{V}}{\sqrt{12}}+\frac{\mathrm{osc}}{\pi }.$$

一般效率估计的机制 ​

由 $J_i=-\int f{\varphi }_i$ 和 $\Vert {\varphi }_i{\Vert }_K^2=h_K/3$，对两侧分别使用 Cauchy–Schwarz，再用 $(a+b{)}^2\le 2a^2+2b^2$，得到

$$w_i{J}_i^2\le \frac{2}{3}(h_i^2\Vert f{\Vert }_{K_i}^2+h_{i+1}^2\Vert f{\Vert }_{K_{i+1}}^2).$$

求和后 $\eta \le \sqrt{7/3}V$。另一方面，上面的误差分解给出 $\overline{V}\le \sqrt{12}E+(\sqrt{12}/\pi )\mathrm{osc}$。结合 $V\le \overline{V}+\mathrm{osc}$，有

$$\eta \le \sqrt{\frac{7}{3}}[\sqrt{12}E+(1+\frac{\sqrt{12}}{\pi })\mathrm{osc}].$$

常数并不尖锐，但不依赖相邻单元长度比。这是含数据振荡的效率，与此前的可靠性方向相反：估计量的偏大受真实误差和未分辨荷载共同限制。[1]

直觉

可靠性不声称 $\eta$ 等于 $E$，也不声称每个 ${\eta }_K$ 与 $E_K$ 成固定比例。相邻单元会共享跳跃贡献，所以真实误差为零的单元也可能收到正指标。

这不是说振荡就是误差，而是说用单元平均值代替原荷载时，遗漏的影响需要另外控制。振荡小，常荷载公式就有解释力；振荡大，单靠平均荷载无法可靠评估误差。

例子与边界

高频荷载揭示的边界 ​

只用一个单元 $(0,1)$，零边界一次元空间是 $V_h=\{0\}$。取 $f_m(x)=\sin \left(2m\pi x\right)$，则

$$u_m(x)=\frac{\sin \left(2m\pi x\right)}{(2m\pi {)}^2},E_m=\frac{1}{2m\pi \sqrt{2}},{\eta }_m=V_m=\frac{1}{\sqrt{2}}.$$

于是 ${\eta }_m/E_m=2m\pi \to \mathrm{\infty }$：本估计量不可能对任意荷载满足一个与 $m$ 无关的无振荡效率界。椭圆方程把高频荷载的响应压小，荷载的 $L^2$ 范数却没有随频率下降。

此时 $\overline{f}=0$、$\mathrm{osc}=1/\sqrt{2}$。若先把荷载替换成单元平均值，再把振荡项删掉，甚至会报告零估计量，而真实误差严格为正。这个反例解释了为什么效率定理中的振荡项具有实际内容。

推论与应用

从误差界走向计算决策 ​

$\eta /\pi \le \epsilon$ 在本条精确积分、精确求解的模型中足以保证 $E\le \epsilon$。若采用近似系数 $\hat{U}$，此前的节点插值恒等式一般不再成立，不能不加修改地套用证明。令 $r=b-A\hat{U}$，则到精确离散解的代数误差满足

$$E_{\mathrm{alg}}^2=r^{\mathsf{T}}{A}^{-1}r.$$

Galerkin 正交还给出总能量误差平方等于离散误差平方与代数误差平方之和。裸残差 $\Vert r{\Vert }_2$ 需要通过逆矩阵或可靠的等价界才能转成这个能量误差。

一旦离散误差占主导，继续求解同一矩阵不会改善网格逼近；应该根据 ${\eta }_K^2$ 选择单元。下一条的完整加密任务公理库自适应有限元的估计、标记与加密Adaptive finite element method · AFEM · Dörfler marking在局部荷载的一维 Poisson 问题上完成求解、残差估计、平方指标标记与局部二分，并比较相同误差和相同预算下的网格选择。会展示：初始四个单元中，只二分有荷载的那个单元，就能取得均匀二分八个单元的相同能量误差，同时也检验“指标较大”与“局部真误差较大”之间并非逐单元恒等。

参考资料

[1] Ricardo H. Nochetto, Kunibert G. Siebert and Andreas Veeser, Theory of Adaptive Finite Element Methods: An Introduction, 2009，作者讲义，定理 13 与式 (81)–(83)，第 80 页；定理 14，第 86 页；推论 9 与注 26，第 88 页。分别讨论残差可靠性、局部效率及数据振荡。本条选定的一维权重和显式常数由正文独立证明。

[2] Ivo Babuška and W. C. Rheinboldt, “Error Estimates for Adaptive Finite Element Computations,” SIAM Journal on Numerical Analysis 15(4), 736–754, 1978，原文。第 3 节以局部化残差组织后验误差估计，是这一思路的早期工作；其问题与网格假设不应直接替代本条的具体证明。