数值方法资格考试问题和解答（马里兰大学）/2008 年 1 月

问题 1

考虑系统$Ax=b\,\!$。GMRES 方法从点$\,\!x_{0}$ 开始，并对残差进行归一化$\,\!r_{0}=b-Ax_{0}$，使得$\,\!v_{1}={\frac {r_{0}}{\nu }}$ 的 2 范数为 1。然后它构建正交的 Krylov 基$\,\!V_{k}=(v_{1}\,v_{2}\,\cdots v_{m})$ 满足 $\,\!AV_{k}=V_{k+1}H_{k}$ 其中$\,\!H_{k}$ 是一个$\,\!(k+1)\times k$ 上 Hessenberg 矩阵。然后人们寻找$\displaystyle x$ 的近似形式 $\,\!x(c)=x_{0}+V_{k}c$ 选择$\,\!c_{k}$ 使$\,\!\|r(c)\|=\|b-Ax(c)\|$ 最小化，其中$\,\!\|\cdot \|$ 是通常的欧几里得范数。

问题 1a

证明$\,\!c_{k}$ 使$\|\nu e_{1}-H_{k}c\|$ 最小化。

解答 1a

我们需要证明$\,\!\|b-Ax(c)\|=\|\nu e_{1}-H_{k}c\|$

${\begin{aligned}\|b-Ax(c)\|&=\|b-A(x_{0}+V_{k}c)\|\\&=\|b-Ax_{0}-AV_{k}c\|\\&=\|r_{0}-AV_{k}c\|\\&=\|r_{0}-V_{k+1}H_{k}c\|\\&=\|\nu v_{1}-V_{k+1}H_{k}c\|\\&=\|V_{k+1}\underbrace {(\nu e_{1}-H_{k}c)} _{h_{c}}\|\\&=(V_{k+1}h_{c},V_{k+1}h_{c})^{\frac {1}{2}}\\&=((V_{k+1}h_{c})^{T}V_{k+1}h_{c})^{\frac {1}{2}}\\&=(h_{c}^{T}V_{k+1}^{T}V_{k+1}h_{c})^{\frac {1}{2}}\\&=(h_{c}^{T}h_{c})^{\frac {1}{2}}\\&=\|h_{c}\|\\&=\|\nu e_{1}-H_{k}c\|\end{aligned}}$

问题 1b

假设我们选择使用正交三角化（QR）方法来求解问题a中的最小二乘问题以求解$c_{k}\!\,$。该方法的浮点运算顺序是多少？请给出原因。

解答 1b

我们希望将$H_{k}\!\,$，一个$(k+1)\times k\!\,$ 上 Hessenberg 矩阵，转换为 QR 形式。

成本大约为$4k^{2}+{\frac {1}{2}}k^{2}=4.5k^{2}$，来自 Given's 旋转和回代的成本。

Given's 旋转成本

我们需要$k\!\,$ 个 Given's 旋转相乘以将每个$k\!\,$ 个次对角线项置零，从而将$H_{k}\!\,$ 转换为上三角形式$R\!\,$，

每次连续的 Given's 旋转乘以一个上 Hessenberg 矩阵需要减少四个乘法，因为每个先前的次对角线项已被 Given's 旋转乘法置零。

因此，$k\!\,$ 个 Given's 旋转相乘的成本为

$4k+4(k-1)+\ldots +4\cdot 1<4k\cdot k=4k^{2}$

回代成本

$R\!\,$ 是一个$(k+1)\times k\!\,$ 上三角矩阵，最后一行全为零。因此，我们需要回代$k\!\,$ 行上三角。

$1+2+\ldots +k={\frac {k(k+1)}{2}}={\frac {k^{2}}{2}}+{\frac {k}{2}}$

问题 2

我们希望近似积分$I(f)\equiv \int _{a}^{b}f(x)dx$

问题 2a

说明复合梯形法则$Q_{T,n}\!\,$ 用于逼近$I(f)\!\,$ 在宽度为$h={\frac {b-a}{n}}\!\,$ 的均匀划分上，并给出误差$I(f)-Q_{T,n}\!\,$ 的公式，该公式适合外推。

解 2a

复合梯形法则为

$Q_{T,n}(f)={\frac {h}{2}}(f(a)+f(b))+h\sum _{k=1}^{n-1}f(x_{k})\!\,$

误差为

$I(f)-Q_{T,n}(f)=-{\frac {(b-a)f^{\prime \prime }(\xi )}{12}}h^{2}$

其中 $\xi \in [a,b]\!\,$.

复合梯形误差推导

局部误差，即在一个区间上的误差，为

$-{\frac {h^{3}}{12}}f^{\prime \prime }(\eta )\!\,$.

观察到

$n\min _{\eta \in [a,b]}f^{\prime \prime }(\eta _{i})\leq \sum _{i=1}^{n}f^{\prime \prime }(\eta _{i})\leq n\max _{\eta \in [a,b]}f^{\prime \prime }(\eta _{i})\!\,$

这意味着

$\min _{\eta \in [a,b]}f^{\prime \prime }(\eta _{i})\leq \sum _{i=1}^{n}{\frac {f^{\prime \prime }(\eta _{i})}{n}}\leq \max _{\eta \in [a,b]}f^{\prime \prime }(\eta _{i})\!\,$

因此，中值定理意味着存在一个 $\xi \in [a,b]\!\,$ 使得

$f^{\prime \prime }(\xi )=\sum _{i=1}^{n}{\frac {f^{\prime \prime }(\eta _{i})}{n}}\!\,$.

将等式两边都乘以 $n\!\,$,

$nf^{\prime \prime }(\xi )=\sum _{i=1}^{n}f^{\prime \prime }(\eta _{i})\!\,$

使用此关系，我们有

${\begin{aligned}I(f)-Q_{T,n}(f)&=\sum _{i=1}^{n}I_{i}(f)-Q_{i}(f)\\&=\sum _{i=1}^{n}\left[-{\frac {h^{3}}{12}}f^{\prime \prime }(\eta _{i})\right]\\&=-{\frac {h^{3}}{12}}f^{\prime \prime }(\xi )n\\&=-{\frac {h^{3}}{12}}f^{\prime \prime }(\xi )({\frac {b-a}{h}})\\&=-{\frac {h^{2}}{12}}f^{\prime \prime }(\xi )(b-a)\end{aligned}}$

局部误差的推导

多项式插值的误差可以通过使用以下定理找到

 Assume $f^{n+1}\!\,$ exists on $[a,b]\!\,$ and $\{x_{0},x_{1},\ldots ,x_{n}|x\in [a,b]\}.\!\,$ $P_{n}\!\,$ interpolates $f\!\,$ at $\{x_{j}\}_{j=0}^{n}\!\,$.  
 Then there is a $\xi _{x}\!\,$ ($\xi \!\,$ is dependent on $x\!\,$) such that 
 
 $f(x)-p_{n}(x)={\frac {(x-x_{0})(x-x_{1})...(x-x_{n})}{(n+1)!}}f^{(n+1)}(\xi _{x})\!\,$
 
 where $\xi _{x}\!\,$ lies in $(m,M)\!\,$,  $m=\min\{x_{0},x_{1},\ldots ,x_{n},X\}\!\,$,  $M=\max\{x_{0},x_{1},\ldots ,x_{n},X\}\!\,$

应用该定理得到：

$f(x)-p_{1}(x)=(x-a)(x-b){\frac {f^{\prime \prime }(\xi _{x})}{2}}$

因此，

${\begin{aligned}E(f)&=\int _{a}^{b}f(x)-p_{1}(x)dx\\&=\int _{a}^{b}\underbrace {\underbrace {(x-a)} _{>0}\underbrace {(x-b)} _{<0}} _{w(x)}{\frac {f^{\prime \prime }(\xi _{x})}{2}}dx\end{aligned}}$

由于 $a\!\,$ 是区间的起点，$x-a\!\,$ 总是正数。相反，由于 $b\!\,$ 是区间的终点，$x-b\!\,$ 总是负数。因此，$w(x)<0\!\,$ 的符号始终一致。因此，根据 积分中值定理，存在一个 $\zeta \in [a,b]\!\,$ 使得

$E(f)={\frac {f^{\prime \prime }(\zeta )}{2}}\int _{a}^{b}(x-a)(x-b)dx\!\,$

请注意，$f^{\prime \prime }(\zeta )$ 是一个常数，与 $x\!\,$ 无关。

对 $\int _{a}^{b}(x-a)(x-b)dx\!\,$ 进行积分，得到：

${\begin{aligned}E(f)&={\frac {f^{\prime \prime }(\zeta )}{2}}\left(-{\frac {(b-a)}{6}}\right)\\&=-{\frac {f^{\prime \prime }(\zeta )(b-a)}{12}}\!\,\end{aligned}}$

问题 2b

使用误差公式推导一个新的求积公式，该公式通过对复合梯形公式进行一次外推得到。这个公式是什么？它的误差如何依赖于 $h\!\,$？你可以假设 $f\!\,$ 是你需要的那样光滑。

解 2b

STOER AND BUESCH 第 162 页

RICHARD HAMMING 的《应用数值分析导论》第 178 页

复合梯形公式在 $n\!\,$ 个点时的误差为

$E_{n}=I(f)-Q_{T,n}={\frac {-(b-a)h^{2}}{12}}+{\mathcal {O}}(h^{3})\!\,$

当有 $2n\!\,$ 个点时，有双倍的区间，但是区间宽度减半。因此，复合梯形公式在 $2n\!\,$ 个点时的误差为

$E_{2n}=I(f)-Q_{T,2n}=2\cdot {\frac {-(b-a)({\frac {h}{2}})^{2}}{12}}+{\mathcal {O}}(h^{3})={\frac {-(b-a)h^{2}}{24}}+{\mathcal {O}}(h^{3})\!\,$

我们可以通过选择 $E_{n}\!\,$ 和 $E_{2n}\!\,$ 的适当线性组合来消除 $h^{2}\!\,$ 项。这将得到一个新的误差规则，该规则的误差为 $h^{3}\!\,$。

${\begin{aligned}E_{n}-2E_{2n}&={\frac {-(b-a)h^{2}}{12}}-2\cdot {\frac {-(b-a)h^{2}}{24}}+{\mathcal {O}}(h^{3})\\&={\mathcal {O}}(h^{3})\end{aligned}}$

将我们对 $E_{n}\!\,$ 和 $E_{2n}\!\,$ 的方程代入左侧，得到

$I(f)-Q_{T,n}-2I(f)-2Q_{T,2n}={\mathcal {O}}(h^{3})\!\,$

$I(f)=2Q_{T,2n}-Q_{T,n}+{\mathcal {O}}(h^{3})\!\,$

如果我们把这个新规则称为 ${\hat {Q}}\!\,$，我们有

${\hat {Q}}=2Q_{T,2n}-Q_{T,n}\!\,$，其误差为 ${\mathcal {O}}(h^{3})\!\,$ 级。

问题 3

考虑用于计算 2 x 2 矩阵 $A\!\,$ 的移位 QR 迭代：从 $A_{0}=A\!\,$ 开始，计算 $A_{k}-\sigma _{k}I=Q_{k}R_{k}\quad \quad \quad A_{k+1}=R_{k}Q_{k}+\sigma _{k}I$ 其中 $\sigma _{k}\!\,$ 是一个标量位移。

问题 3a

如果 $A_{k}={\begin{bmatrix}a_{11}&a_{12}\\a_{21}&a_{22}\end{bmatrix}},$ 指定正交矩阵 $Q_{k}\!\,$，在使用 Givens 旋转进行此步骤时，该矩阵应使用平移矩阵的条目进行描述。

解 3a

设 ${\tilde {A}}_{k}\!\,$ 是 $\sigma _{k}\!\,$ 平移的 $A_{k}\!\,$ 矩阵，即

${\tilde {A}}_{k}=A_{k}-\sigma _{k}I={\begin{bmatrix}a_{11}-\sigma _{k}&a_{12}\\a_{21}&a_{22}-\sigma _{k}\end{bmatrix}},$

$Q_{k}^{T}\!\,$ 是一个 2x2 的正交 Givens 旋转矩阵。 $Q_{k}^{T}\!\,$ 的条目被选择，这样当 $Q_{k}^{T}\!\,$ 乘以 2x2 矩阵 ${\tilde {A}}_{k}\!\,$ 时， $Q_{k}^{T}\!\,$ 将使 ${\tilde {A}}_{k}\!\,$ 的 (2,1) 项为零，并按比例缩放 ${\tilde {A}}_{k}\!\,$ 的剩余三个条目，即

$Q_{k}^{T}{\tilde {A}}_{k}={\begin{bmatrix}*&*\\0&*\end{bmatrix}}=R_{k}$

其中 $*\!\,$ 表示我们不感兴趣的可计算标量值，而 $R_{k}\!\,$ 是我们使用 QR 算法寻求的所需上三角矩阵。

由于 $Q_{k}^{T}\!\,$ 是正交的，上面的等式意味着

${\begin{aligned}Q_{k}^{T}{\tilde {A}}_{k}&=R_{k}\\Q_{k}Q_{k}^{T}{\tilde {A}}_{k}&=Q_{k}R_{k}\\{\tilde {A}}_{k}&=Q_{k}R_{k}\end{aligned}}$

给定旋转 $Q_{k}^{T}\!\,$ 由下式给出

$Q_{k}^{T}={\begin{bmatrix}c&s\\-s&c\end{bmatrix}}$

其中

${\begin{aligned}c&={\frac {a_{11}-\sigma _{k}}{r}}\\s&={\frac {a_{21}}{r}}\\r&={\sqrt {(a_{11}-\sigma _{k})^{2}+a_{21}^{2}}}\end{aligned}}$

对 $Q_{k}^{T}\!\,$ 求转置，得到

$Q_{k}={\begin{bmatrix}c&-s\\s&c\end{bmatrix}}$

问题 3b

假设 $a_{21}=\delta \!\,$，一个小数，并且 $|a_{12}|\leq (a_{11}-a_{22})^{2}\!\,$。证明，通过适当的平移 $\sigma _{k}\!\,$，$A_{k+1}\!\,$ 的 (2,1) 项的大小为 $O(\delta ^{2})\!\,$。这对于 QR 迭代的收敛速度说明了什么？

解 3b

根据假设和（a）部分，

${\begin{aligned}A_{k+1}&=R_{k}Q_{k}+\sigma _{k}I\\&={\begin{bmatrix}*&*\\0&r_{22}\end{bmatrix}}{\begin{bmatrix}c&-s\\s&c\end{bmatrix}}+{\begin{bmatrix}\sigma _{k}&0\\0&\sigma _{k}\end{bmatrix}}\end{aligned}}$

令 $A_{k+1}(2,1)\!\,$ 是 $A_{k}\!\,$ 的 (2,1) 项。利用上面的等式，我们可以通过求 $R_{k}\!\,$ 的第二行和 $Q_{k}\!\,$ 的第一列的内积，并加上 $\sigma _{k}I\!\,$ 的 (2,1) 项来找到 $A_{k+1}(2,1)\!\,$，即

${\begin{aligned}A_{k+1}(2,1)&=(0\cdot c+r_{22}\cdot s)+0\\&=r_{22}s\end{aligned}}$

我们需要找到 $r_{22}\!\,$ 的值，所以我们需要计算 $R_{k}\!\,$.

根据假设和 $Q_{k}\!\,$ 的正交性，我们有

${\begin{aligned}A_{k}-\sigma _{k}I&=Q_{k}R_{k}\\Q_{k}^{T}(A_{k}-\sigma _{k}I)&=Q_{k}^{T}Q_{k}R_{k}\\Q_{k}^{T}(A_{k}-\sigma _{k}I)&=R_{k}\\R_{k}&=Q_{k}^{T}(A_{k}-\sigma _{k}I)\\&={\begin{bmatrix}c&s\\-s&c\end{bmatrix}}{\begin{bmatrix}a_{11}-\sigma _{k}&a_{12}\\\delta &a_{22}-\sigma _{k}\end{bmatrix}}\end{aligned}}$

从 $R_{k}\!\,$ 中，我们可以通过使用内积找到其 (2,2) 元素 $r_{2}2\!\,$。

$r_{22}=-s\cdot a_{12}+c\cdot (a_{22}-\sigma _{k})\!\,$

现在我们有了 $r_{22}\!\,$，我们可以通过使用适当的代入来计算 $A_{k+1}(2,1)\!\,$。

${\begin{aligned}A_{k+1}(2,1)&=r_{22}\cdot s\\&=(-sa_{12}+c(a_{22}-\sigma _{k}))s\\&=-s^{2}a_{12}+cs(a_{22}-\sigma _{k})\\&={\frac {-\delta ^{2}a_{12}}{(a_{11}-\sigma _{k})^{2}+\delta ^{2}}}+{\frac {(a_{11}-\sigma _{k})(a_{22}-\sigma _{k})\delta }{(a_{11}-\sigma _{k})^{2}+\delta ^{2}}}\\&={\frac {-\delta ^{2}a_{12}+\delta (a_{11}-\sigma _{k})(a_{22}-\sigma _{k})}{(a_{11}-\sigma _{k})^{2}+\delta ^{2}}}\\&\approx {\frac {-\delta ^{2}a_{12}+\delta (a_{11}-\sigma _{k})(a_{22}-\sigma _{k})}{(a_{11}-\sigma _{k})^{2}}}\end{aligned}}$

因为 $\delta \!\,$ 是一个很小的数。

令我们的平移 $\sigma _{k}=a_{22}\!\,$。 那么上面的等式将得到：

${\begin{aligned}A_{k+1}(2,1)&\approx {\frac {-\delta ^{2}a_{12}+\delta (a_{11}-\sigma _{k})(a_{22}-\sigma _{k})}{(a_{11}-\sigma _{k})^{2}}}\\&\approx {\frac {-\delta ^{2}a_{12}+\delta (a_{11}-a_{22})(a_{22}-a_{22})}{(a_{11}-a_{22})^{2}}}\\&={\frac {-\delta ^{2}a_{12}}{(a_{11}-a_{22})^{2}}}\end{aligned}}$

因此，

${\begin{aligned}|A_{k+1}(2,1)|&\approx \left|{\frac {\delta ^{2}a_{12}}{(a_{11}-a_{22})^{2}}}\right|\\&\leq |\delta ^{2}|\end{aligned}}$

因为 $|a_{12}|\leq (a_{11}-a_{22})^{2}\!\,$

因此，我们已经证明 $A_{k+1}(2,1)\!\,$ 是 $\delta ^{2}\!\,$ 阶。

如果 $\delta \!\,$ 很小，则 QR 收敛速度将是二次的。