自适应有限元的估计、标记与加密

均匀加密把计算资源平均分给整个区域；自适应有限元则先观察当前解在哪里还不够准确，再把新自由度放到这些位置。它的关键不是让网格看起来不均匀，而是让“误差证据—单元选择—重新求解”形成可核查的闭环。

本条完成一个具体任务：从四单元网格出发，计算每个残差指标，选出需要二分的单元，重解系统，并用已知精确解检查误差变化。最后用相同误差、相同自由度两种比较说明加密位置的作用。所有积分和小型线性系统都按精确算术处理。

形式陈述

估计量和标记规则

沿用一维残差后验估计的零端点 Poisson 模型、连续一次元及指标

\eta_K^2=h_K^2\|f\|_K^2+ \frac12\sum_{x_i\in\partial K\cap(0,1)}\min(h_i,h_{i+1})J_i^2, \qquad \eta^2=\sum_K\eta_K^2.

一次循环由四项操作组成：在当前网格上求解离散方程；计算各 $\eta_K^2$ ；选择标记集合 $\mathcal M$ ；二分标记单元并在新空间中重新求解。每个新循环都重新计算指标，不能永久沿用第一轮的排序。

这里采用平方形式的 Dörfler 标记：给定 $0<\theta\le1$ ，要求

\sum_{K\in\mathcal M}\eta_K^2\ge\theta\eta^2.

为使任务完全确定，先按 $\eta_K^2$ 从大到小排序，相等时从左到右；取达到门槛的最短前缀。它具有最小的标记单元数，因为任意同样数量的单元都不会超过最大若干项的总和。若 $\eta=0$ ，直接停止。

这里 $\theta$ 表示要覆盖的平方指标总量比例。有些文献写成 $\eta(\mathcal M)\ge\vartheta\eta$ ，其中 $\eta(\mathcal M)^2=\sum_{K\in\mathcal M}\eta_K^2$ ；两者对应 $\theta=\vartheta^2$ ，不能把参数数值直接混用。[1]

先明确停止证书

对精确离散解，已有 $E\le\eta/\pi$ ，所以可以用 $\eta/\pi\le\varepsilon$ 作为本模型的充分停止条件。指标下降是计算进展，达到具有可靠性常数的门槛才构成误差保证。实际迭代求解的额外误差需要另行控制，不能用一个未经解释的矩阵残差阈值替代。

直觉

局部指标把连续方程尚未满足的部分分配到单元上；标记规则再把这些分散证据转成有限的加密选择。覆盖一半平方指标总量，并不要求标记一半单元：如果误差证据集中在少量位置，只增加少量节点就可能显著改善逼近。重新求解后，斜率和跳跃都会改变，因此下一轮应依据新指标重新作决定。

一维二分只需在被选区间中加入中点，相邻区间仍然形成协调划分。高维加密还涉及邻居配合、悬挂节点与形状控制，这些操作不由本条的一维循环自动涵盖。

例子与边界

四单元网格上的求解与估计

局部荷载让误差集中在左端

取

f(x)=\begin{cases}4,&0<x<1/4,\\0,&1/4<x<1,\end{cases} \qquad u(0)=u(1)=0.

精确解为

u(x)=\begin{cases} 7x/8-2x^2,&0\le x\le1/4,\\ (1-x)/8,&1/4\le x\le1. \end{cases}

两段在 $x=1/4$ 的函数值都是 $3/32$ ，导数都是 $-1/8$ ，因此解及通量连续。左端是抛物线，右侧已经是一次函数。这意味着一次元需要改善的弯曲只存在于左端；右侧即使单元很长也能精确表达这一段解。

初始节点为 $0,1/4,1/2,3/4,1$ 。荷载断点与节点对齐，因此每个单元的数据振荡为零，局部真误差可以使用 $E_K^2=c_K^2h_K^3/12$ 精确核验。

组装、求解并计算跳跃

消去边界值后，系统为

A=4\begin{pmatrix}2&-1&0\\-1&2&-1\\0&-1&2\end{pmatrix},\qquad b=\begin{pmatrix}1/2\\0\\0\end{pmatrix},\qquad U=\begin{pmatrix}3/32\\1/16\\1/32\end{pmatrix}.

第一个荷载分量来自 $4\int_0^{1/4}\phi_1\,dx=1/2$ 。代入矩阵可逐行验证 $AU=b$ 。从节点值计算的四个斜率和三个跳跃是

s=(3/8,-1/8,-1/8,-1/8),\qquad J=(-1/2,0,0).

第一个单元的体积项为

h_1^2\|f\|_{K_1}^2=(1/4)^2\cdot16\cdot(1/4)=1/4.

唯一非零跳跃的加权平方是 $(1/4)(1/2)^2=1/16$ ，左右各收到 $1/32$ 。所以

(\eta_{K_1}^2,\ldots,\eta_{K_4}^2)=(9/32,1/32,0,0),\qquad \eta^2=5/16.

真实误差只在第一单元：

(E_{K_1}^2,\ldots,E_{K_4}^2)=(1/48,0,0,0),\qquad E=1/\sqrt{48}\approx0.144338.

第二单元的真误差为零，但因共享左端跳跃而有正指标。后验指标为选择提供证据，并不逐单元复制真实误差。

标记、二分与重新求解

第一轮只需选一个单元

设 $\theta=1/2$ ，需覆盖的平方指标总量是 $5/32$ 。最大项 $9/32$ 单独就达到要求，因此 $\mathcal M=\{K_1\}$ 。二分后节点变为

0,\quad1/8,\quad1/4,\quad1/2,\quad3/4,\quad1.

新系统与解分别是

A_+=\begin{pmatrix} 16&-8&0&0\\ -8&12&-4&0\\ 0&-4&8&-4\\ 0&0&-4&8 \end{pmatrix},\qquad b_+=\begin{pmatrix}1/2\\1/4\\0\\0\end{pmatrix},\qquad U_+=\begin{pmatrix}5/64\\3/32\\1/16\\1/32\end{pmatrix}.

这次有四个内部自由度。新增节点不是简单取旧折线的中点值；旧插值在那里取 $3/64$ ，重新求解后得到 $5/64$ 。空间增加之后，求解才能利用新增方向修正原来的误差。

新指标与新真误差

新斜率与跳跃为

s_+=(5/8,1/8,-1/8,-1/8,-1/8),\qquad J_+=(-1/2,-1/4,0,0).

前两个单元的体积项各为 $16(1/8)^3=1/32$ 。两个非零节点项分别为 $1/32$ 和 $1/128$ ，仍各半分给左右单元，因此

(\eta_K^2)_+=\left(\frac3{64},\frac{13}{256},\frac1{256},0,0\right),\qquad \eta_+^2=\frac{13}{128}.

两个带荷载的单元各有 $E_K^2=16(1/8)^3/12=1/384$ ，总计

E_+^2=1/192,\qquad E_+=1/\sqrt{192}\approx0.072169.

能量误差恰好减半，而 $\eta_+/\eta=\sqrt{13/40}\approx0.5701$ 。估计量与真误差并没有同比下降，可靠性和效率也没有要求它们同比。

如果再运行一轮，同样的标记门槛是 $\eta_+^2/2=13/256$ 。现在最大项是第二单元的 $13/256$ ，恰好单独达到门槛，因此下一步应二分 $(1/8,1/4)$ ，不是机械地继续二分最左单元。二分后总误差平方为 $1/384+1/1536=5/1536$ ；这由单元误差的四分之一规律得到。

初始四单元网格的误差集中在左端；二分左端后四个内部自由度达到均匀八单元、七个内部自由度的相同误差，二分右端则误差不变。

加密位置决定新增自由度是否能改善本例的误差。蓝色粗段标出荷载所在区域，新增中点以空心圆表示； $E^2$ 是能量误差平方，DOF 是内部自由度数。左端二分与均匀二分具有相同误差，但分别只需 4 与 7 个内部自由度。

相同误差与相同预算的比较

将初始全部四个单元都二分，会得到八个单元。其前两个单元与局部加密网格相同；剩余区间的精确解本来就是直线，多加节点不会降低已经为零的误差。因此均匀二分也有 $E^2=1/192$ ，并且本例中全局估计量也恰好相同。

反过来，若只二分最后一个单元 $(3/4,1)$ ，同样使用五个单元、四个内部自由度，却没有改善左端的逼近。新增节点值为 $u(7/8)=1/64$ ，其相邻斜率仍为 $-1/8$ ，所有非零残差贡献仍在原位置。

网格选择	单元数	内部自由度	$\eta^2$	$E^2$
初始四等分	4	3	$5/16$	$1/48$
只二分左端单元	5	4	$13/128$	$1/192$
所有单元均匀二分	8	7	$13/128$	$1/192$
只二分右端单元	5	4	$5/16$	$1/48$

这张表的结论是具体的：对该荷载，第一次指标标记用较少自由度取得了均匀加密的同一误差；相同预算放错位置则没有改善。它不证明任意荷载、任意标记和任意网格族上的最优复杂度。一般自适应收敛理论还需要分析估计量缩减、网格性质、振荡及相邻离散解之间的关系。[1][2]

代数求解误差不能混进加密判断

能量正交给出的精确分解

设 $U$ 是某张网格的精确离散系数， $\widehat U$ 是近似求解结果， $r=b-A\widehat U$ 。令 $d=U-\widehat U=A^{-1}r$ ，则对应有限元函数之间的能量距离是

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

由于 $u_h-\widehat u_h\in V_h$ ，Galerkin 正交消去交叉项，得到

\|u'-\widehat u_h'\|_2^2 =\|u'-u_h'\|_2^2+\|u_h'-\widehat u_h'\|_2^2 =E_{\mathrm{disc}}^2+E_{\mathrm{alg}}^2.

因此两种工作解决不同部分：改善线性求解消除 $E_{\mathrm{alg}}$ ，改变网格才可能改善 $E_{\mathrm{disc}}$ 。用 $\|r\|_2$ 代替 $r^{\mathsf T}A^{-1}r$ 会丢掉矩阵对残差方向和尺度的作用。

在初始网格上核验一次

故意取 $\widehat U=U+(1/32,0,0)^{\mathsf T}$ ，则

r=(-1/4,1/8,0)^{\mathsf T},\qquad E_{\mathrm{alg}}^2=(1/32)^2A_{11}=1/128.

此时总能量误差平方是 $1/48+1/128=11/384$ 。把系统精确解完后， $r$ 和代数误差都变为零，但离散误差仍为 $1/48$ 。这个计算也解释了为什么本条先明确“精确求解”，再使用残差指标的证明；对于近似离散解，要额外估计代数误差，不能直接沿用节点精确性。

可复算任务的完成标准

可以独立实现这四张网格的组装与求解。每个单元累加 $h_K^{-1}\begin{pmatrix}1&-1\\-1&1\end{pmatrix}$ ，常荷载 $c_K$ 的局部右端为 $(c_Kh_K/2)(1,1)^{\mathsf T}$ ；消去两端自由度后求解。用所得节点值计算斜率、跳跃、局部指标，再与表中分数比较。

真误差必须在单元内部积分 $\int_K(u'-s_K)^2$ ，不能只比较节点值，因为本模型在节点处恰好精确。完成后应能解释四件事：为什么第一轮只标记左端；为什么下一轮改标第二单元；为什么均匀加密浪费了本例右侧的新增自由度；为什么矩阵残差为零仍有正的 PDE 误差。这些解释与分数计算共同构成一次完整的自适应任务。

推论与应用

一次循环实际需要多少计算

设当前网格有 $N$ 个单元。在各单元所需的荷载积分与 $\|f\|_K^2$ 已提供、每次标量运算与比较按常数成本计的模型下，组装、计算斜率与全部局部指标、一维二分及节点重编号各需 $O(N)$ 操作。消去端点后，刚度矩阵是对称正定三对角矩阵，三对角 $LDL^{\mathsf T}$ 分解与回代需 $O(N)$ 标量运算，存储需 $O(N)$ ；二分后单元数至多为 $2N$ ，重解仍具有同一线性量级。本条指定的完整排序需 $O(N\log N)$ 比较，之后扫描前缀需 $O(N)$ ，所以这一实现的单轮成本由排序主导，辅助存储可保持 $O(N)$ 。

这些是算术操作计数，不把精确有理数的位数增长算作常数。对任意 $f$ ，取得精确单元积分也不是免费步骤；若改用数值积分，其计算量及引入的误差须另外计算。经过多轮、各轮网格大小为 $N_0,\ldots,N_L$ 时，应累加各轮成本，而不能仅用最后一张网格的大小代替总工作量。前面的同误差比较只核验了这个算例的自由度收益，没有证明达到给定误差所需的总计算量最优。

参考资料

[1] Ricardo H. Nochetto, Kunibert G. Siebert and Andreas Veeser, Theory of Adaptive Finite Element Methods: An Introduction, 2009，作者讲义。第 93 页的 Mark 条件使用估计量范数形式的参数；§8.1 与引理 18，第 109–110 页，讨论嵌套空间、能量正交及自适应分析。本文的平方参数已明确换算，所列数值由当前一维模型直接计算。

[2] Willy Dörfler, “A Convergent Adaptive Algorithm for Poisson’s Equation,” SIAM Journal on Numerical Analysis 33(3), 1106–1124, 1996，论文页面。这是 Dörfler 标记的原始工作；本条不把该论文特定算法的收敛结论直接赋予上述一维实现。

形式陈述 ​

估计量和标记规则 ​

先明确停止证书 ​

四单元网格上的求解与估计 ​

局部荷载让误差集中在左端 ​

组装、求解并计算跳跃 ​

标记、二分与重新求解 ​

第一轮只需选一个单元 ​

新指标与新真误差 ​

相同误差与相同预算的比较 ​

代数求解误差不能混进加密判断 ​

能量正交给出的精确分解 ​

在初始网格上核验一次 ​

可复算任务的完成标准 ​

一次循环实际需要多少计算 ​

反向链接