数值方法资格考试试题及解答 (马里兰大学)/Jan05 667

问题 4

设 $f:R^{n}\rightarrow R\!\,$ 是一个非线性光滑函数。为了确定 $f\!\,$ 的（局部）最小值，可以使用如下形式的下降方法

x_{k+1}=x_{k}+\alpha _{k}d_{k}\quad (k\geq 0)\quad \quad (1)\!\,

其中 $0\leq \alpha _{k}\leq 1\!\,$ 是通过回溯获得的适当参数，而 $d_{k}$ 是一个下降方向，即它满足

f(x_{k+1})<f(x_{k})\quad \quad (2)\!\,

问题 4a

写出最速下降法（或梯度法），并证明存在 $0<\alpha _{k}\leq 1\!\,$ 使得由此产生的方法满足 (2)

解答 4a

最速下降法

$x_{k+1}=x_{k}-\alpha _{k}\nabla f(x_{k})\!\,$

选择 $\alpha _{k}\in (0,1]\!\,$ 使得 $\alpha _{k}\!\,$ 最小化 $f(x_{k}-\alpha _{k}\nabla f(x_{k}))\!\,$ ，即

$\alpha _{k}:\min _{0<\alpha _{k}\leq 1}f(x_{k}-\alpha _{k}\nabla f(x_{k}))\!\,$

方向导数为负

为了满足 (2)，方向导数应该为负，即

$\nabla _{d_{k}}f(x_{k})=\lim _{\alpha _{k}\rightarrow 0}{\frac {f(x_{k}+\alpha _{k}d_{k})-f(x_{k})}{\alpha _{k}}}<0\!\,$

这意味着

$\underbrace {f(x_{k}+\alpha _{k}d_{k})} _{f(x_{k+1})}<f(x_{k})\!\,$

因为 $\alpha _{k}\in (0,1]\!\,$

问题 4b

写出牛顿法，并检验是否存在 $\alpha _{k}\!\,$ 使得 (2) 成立。确定 $f(x)\!\,$ 的海森矩阵 $H(f(x))\!\,$ 上的条件，这些条件保证了 $\alpha _{k}\!\,$ 的存在。

解 4b

牛顿法

$x_{k+1}=x_{k}-H^{-1}(f(x_{k}))\nabla f(x_{k})\!\,$

方向导数为负

为了下降，我们需要方向导数为负，即

$\nabla _{\alpha _{k}d_{k}}f(x_{k})=\nabla f(x_{k})\cdot \alpha _{k}d_{k}=-\nabla f(x_{k})\cdot H^{-1}(f(x_{k}))\nabla f(x_{k})<0\!\,$

这意味着

$\nabla f(x_{k})^{T}H^{-1}(f(x_{k}))\nabla f(x_{k})>0\!\,$

因此 $H^{-1}(f(x_{k}))\!\,$ 是正定的，因此 $H(f(x_{k}))\!\,$ 是正定的。

问题 4c

如果我们用矩阵 $H(f(x))+\gamma _{k}I\!\,$ 替换 Hessian，其中 $\gamma _{k}\geq 0\!\,$ 且 $I$ 是单位矩阵，我们得到一个拟牛顿法。求 $\gamma _{k}\!\,$ 的条件，它会导致 (2)。

解 4c

我们现在需要 ${\tilde {H}}(f(x_{k}))=H(f(x_{k}))+\gamma _{k}I\!\,$ 是正定的。

令 $\lambda _{i}\!\,$ ， $i=1,2,\ldots ,n\!\,$ 是 $H(f(x_{k}))\!\,$ 的特征值。

然后 $\lambda _{i}+\gamma _{k}\cdot 1\!\,$ ， $i=1,2,\ldots ,n\!\,$ 是 ${\tilde {H}}(f(x_{k}))\!\,$ 的特征值。

由于我们希望 ${\tilde {H}}(f(x_{k}))\!\,$ 为正定，等价于对 $i=1,2,\ldots ,n\!\,$

$\lambda _{i}+\gamma _{k}>0\!\,$

即

$\gamma _{k}>-\lambda _{i}\quad i=1,2,\ldots ,n\!\,$

问题 5

考虑以下在 $R$ 中的非线性自治初值问题，其中 $f\in C^{2}$ 。

y'=f(y)\quad \quad y(0)=y_{0}

问题 5a

将 ODE 写成积分形式，并使用中点求积规则推导出具有均匀时间步长 $h={\frac {T}{N}}$ 的中点法。

y_{n+1}=y_{n}+hf({\frac {y_{n+1}+y_{n}}{2}})

解答 5a

对于 $x\in [x_{n},x_{n+1}]\!\,$ ，

$y(x_{n+1})-y(x_{n})=\int _{x_{n}}^{x_{n+1}}f(y(x))dx\approx h\cdot f\left({\frac {y_{n+1}+y_{n}}{2}}\right)\!\,$

问题 5b

定义截断误差 $T_{n}\!\,$ . 假设 $T_{n}=O(h^{3})\!\,$ , 证明误差 $|y(t_{N})-y_{N}|\!\,$ 的估计值。该方法的阶数是多少？（提示：利用 $f\!\,$ 是 Lipschitz 连续的）。

解 5b

${\begin{aligned}|e_{n+1}|&=|y(t_{n+1})-y_{n+1}|\\&\leq |y(t_{n})-y_{n}|+h\underbrace {\left|f\left({\frac {y(t_{n+1})-y(t_{n})}{2}}\right)-f\left({\frac {y_{n+1}-y_{n}}{2}}\right)\right|} _{\leq {\frac {\gamma }{2}}|y(t_{n+1})-y(t_{n})-y_{n+1}+y_{n}|}+{\mathcal {O}}(h^{3})\\&\leq |e_{n}|+{\frac {h\gamma }{2}}\left(|e_{n+1}|+|e_{n}|\right)+{\mathcal {O}}(h^{3})\end{aligned}}$

其中 $\gamma \!\,$ 是 $f\!\,$ 的 Lipschitz 常数。重新排列项，我们得到

$|e_{n+1}|\leq \underbrace {\frac {1+{\frac {h\gamma }{2}}}{1-{\frac {h\gamma }{2}}}} _{C}|e_{n}|+{\mathcal {O}}(h^{3})\!\,$

特别是， $|e_{1}|\leq C\underbrace {|e_{0}|} _{0}+{\mathcal {O}}(h^{3})\!\,$

那么 $e_{N}\!\,$ 由下式给出

${\begin{aligned}e_{N}&=\left(C^{N-1}+C^{N-2}+\cdots +C+1\right){\mathcal {O}}(h^{3})\\&={\frac {C^{N}-1}{C-1}}{\mathcal {O}}(h^{3})\\&\leq K{\mathcal {O}}(h^{2})\end{aligned}}\!\,$

问题 5c

证明 $T_{n}=O(h^{3})\!\,$ 对于此方法。(提示：首先展开 $y(t_{n+1})\!\,$ 和 $y(t_{n})\!\,$ 在 $\tau _{n}=t_{n}+{\frac {h}{2}}\!\,$ 附近，然后展开 $f({\frac {1}{2}}(y(t_{n+1})+y(t_{n}))).\!\,$ 也展开 $y'(t)\!\,$ 在 $\tau _{n}\!\,$ 附近。）

解 5c

中点法可以改写如下

$y(t_{n+1})=y(t_{n})+{\frac {h}{2}}y'(t_{n+1})+{\frac {h}{2}}y'(t_{n})+T_{n}\!\,$

这意味着

$y(t_{n+1})-y(t_{n})-{\frac {h}{2}}y'(t_{n+1})-{\frac {h}{2}}y'(t_{n})=T_{n}\!\,$

将每项在 $\tau _{n}=t_{n}+{\frac {h}{2}}\!\,$ 附近展开，得到 $T_{n}\!\,$ .

${\begin{array}{|c|c|c|c|c|c|}{\mbox{Order }}h&y(t_{n+1})&-y(t_{n})&-{\frac {h}{2}}y'(t_{n+1})&-{\frac {h}{2}}y'(t_{n})&\Sigma \\\hline &&&&&\\0&y(t_{n}+{\frac {h}{2}})&-y(t_{n}+{\frac {h}{2}})&---&---&0\\&&&&&\\1&y'(t_{n}+{\frac {h}{2}}){\frac {h}{2}}&-y'(t_{n}+{\frac {h}{2}})(-{\frac {h}{2}})&-{\frac {h}{2}}y'(t_{n}+{\frac {h}{2}})&-{\frac {h}{2}}y'(t_{n}+{\frac {h}{2}})&0\\&&&&&\\2&{\frac {y''(t_{n}+{\frac {h}{2}})}{2}}{\frac {h^{2}}{4}}&-{\frac {y''(t_{n}+{\frac {h}{2}}}{2}}{\frac {h^{2}}{4}}&-{\frac {h}{2}}y''(t_{n}+{\frac {h}{2}}){\frac {h}{2}}&-{\frac {h}{2}}y''(t_{n}+{\frac {h}{2}})(-{\frac {h}{2}})&0\\&&&&&\\3&O(h^{3})&O(h^{3})&O(h^{3})&O(h^{3})&O(h^{3})\\\hline \end{array}}$

问题 6

考虑以下边界为两点的边值问题，在 $(0,1)\!\,$ 内，且 $\epsilon <<1\leq b:\!\,$

-\epsilon u_{xx}+bu_{x}=1,\quad u(0)=u(1)=0\!\,

问题 6a

写出有限元方法，其中分段线性元素在 $(0,1)\!\,$ 上形成一个均匀网格，网格尺寸为 $h\!\,$ 。如果 $U={U_{i}}_{i=0}^{I+1}\!\,$ 是有限元解的节点值向量，找到（刚度）矩阵 $A\!\,$ 和右手边 $F\!\,$ ，使得 $AU=F\!\,$ 。 $A\!\,$ 是否对称？ $A\!\,$ 是否正定？

解答 6a

弱变分公式

用测试函数 $v\in H_{0}^{1}[0,1]\!\,$ 乘以给定方程，并从 0 到 1 进行积分，得到等效的弱变分公式

找到 $u\in H_{0}^{1}[0,1]\!\,$ ，使得对于所有 $v\in H_{0}^{1}[0,1]\!\,$ ，以下成立

\underbrace {\epsilon \int _{0}^{1}u'v'dx+b\int _{0}^{1}u'vdx} _{a(u,v)}=\underbrace {\int _{0}^{1}vdx} _{f(v)}\!\,

离散变分公式

令 $V_{h}=\{v\in H_{0}^{1}[0,1]:\left.v_{h}\right|_{[x_{i-1},x_{i}]}{\mbox{ linear }}\!\,$ 当 $i=1,2,\ldots ,n\}\!\,$

那么我们有离散变分公式，它是弱变分公式的近似。

寻找 $u_{h}\in V_{h}\!\,$ 使得对于所有的 $v_{h}\in V_{h}\!\,$

a(u_{h},v_{h})=f(v_{h})\!\,

定义 V_h 的基底

令 $\{\phi _{i}\}_{i=1}^{n}\!\,$ 是线性“帽”函数，它定义了 $V_{h}\!\,$ 的基底。

$\phi _{i}(x)={\begin{cases}{\frac {x-x_{i-1}}{h}}&{\mbox{ if }}x\in [x_{i-1},x_{i}]\\{\frac {x_{i+1}-x}{h}}&{\mbox{ if }}x\in [x_{i},x_{i+1}]\\0&{\mbox{ otherwise }}\end{cases}}\!\,$

然后计算得到以下结果：（画图）

$\int _{0}^{1}\phi _{i}(x)dx=h\!\,$

$\phi _{i}'(x)={\begin{cases}{\frac {1}{h}}&{\mbox{ if }}x\in [x_{i-1},x_{i}]\\-{\frac {1}{h}}&{\mbox{ if }}x\in [x_{i},x_{i+1}]\\0&{\mbox{ otherwise }}\end{cases}}\!\,$

$\int _{0}^{1}\phi _{i}'(x)\phi _{j}'(x)dx={\begin{cases}{\frac {2}{h}}&{\mbox{ if }}i=j\\-{\frac {1}{h}}&{\mbox{ if }}|i-j|=1\\0&{\mbox{ otherwise }}\end{cases}}\!\,$

$\int _{0}^{1}\phi _{i}(x)\phi _{j}'(x)dx={\begin{cases}0&{\mbox{ if }}i=j\\{\frac {1}{2}}&{\mbox{ if }}i-j=1\\-{\frac {1}{2}}&{\mbox{ if }}i-j=-1\\0&{\mbox{ otherwise }}\end{cases}}\!\,$

矩阵形式的离散变分公式

由于 $\{\phi _{i}\}_{i=1}^{n}\!\,$ 是 $V_{h}\!\,$ 的一个基，

u_{h}=\sum _{i=1}^{n}u_{hi}\phi _{i}\!\,

此外，离散变分公式也可以表示为

\epsilon \sum _{i=1}^{n}u_{hi}\int _{0}^{1}\phi _{i}'(x)\phi _{j}'(x)dx+b\sum _{i=1}^{n}u_{hi}\int _{0}^{1}\phi _{i}'(x)\phi _{j}(x)dx=\int _{0}^{1}\phi _{j}(x)dx{\mbox{ for }}j=1,2,\ldots ,n\!\,

其矩阵形式为

$\underbrace {\begin{bmatrix}2{\frac {\epsilon }{h}}&-{\frac {\epsilon }{h}}+{\frac {b}{2}}&0&\cdots &\cdots &0\\-{\frac {\epsilon }{h}}-{\frac {b}{2}}&2{\frac {\epsilon }{h}}&-{\frac {\epsilon }{h}}+{\frac {b}{2}}&0&\cdots &0\\0&\ddots &\ddots &\ddots &\cdots &0\\0&\cdots &\cdots &0&-{\frac {\epsilon }{h}}-{\frac {b}{2}}&2{\frac {\epsilon }{h}}\end{bmatrix}} _{A}{\begin{bmatrix}u_{1}\\u_{2}\\\vdots \\u_{n}\end{bmatrix}}={\begin{bmatrix}h\\h\\\vdots \\h\end{bmatrix}}\!\,$

$A\!\,$ 不是对称矩阵。 $A\!\,$ 如果满足以下条件，则为正定矩阵

{\frac {\epsilon }{h}}>{\frac {b}{2}}\!\,

问题 6b

找出三个参数 $h,\epsilon ,b\!\,$ 之间的关系，使得 $A\!\,$ 为 $M-\!\,$ 矩阵，即当 $i\neq j\!\,$ 时，有 $a_{ii}>0,a_{ij}\leq 0\!\,$ ，以及 $a_{ii}\geq -\sum _{i\neq j}a_{ij}\!\,$

解 6b

第一行、第二行和最后一行都得出了关于 $A\!\,$ 为 $M\!\,$ - 矩阵的相同不等式。

{\frac {\epsilon }{h}}>{\frac {b}{2}}\!\,

问题 6c

考虑常微分方程的迎风修改

-(\epsilon +{\frac {b}{2}}h)u_{xx}+bu_{x}=1\!\,

证明由此产生的矩阵 $A\!\,$ 是一个 M 矩阵，不受 $h,\epsilon \!\,$ 和 $b\!\,$ 的限制。

解 6c

将 $\epsilon +{\frac {b}{2}}h\!\,$ 代入 $\epsilon \!\,$ 得到以下矩阵 $A\!\,$

${\begin{bmatrix}2{\frac {\epsilon }{h}}+b&-{\frac {\epsilon }{h}}&0&\cdots &\cdots &0\\-{\frac {\epsilon }{h}}-b&2{\frac {\epsilon }{h}}+b&-{\frac {\epsilon }{h}}&0&\cdots &0\\0&\ddots &\ddots &\ddots &\cdots &0\\0&\cdots &\cdots &0&-{\frac {\epsilon }{h}}-b&2{\frac {\epsilon }{h}}+b\end{bmatrix}}\!\,$