清华大学贾忠孝高值第三次作业参考答案

A
= 0 最小, 即近似解与精确解
的误差在 A 范数意义下最小, 算法一步收敛(但是实际中这个 v 不能精确找
解. 设 p = r + αAr, (r, Ap) = 0 等 价 于 (Ar, p) = 0. 解 得 α =
−(Ar, r)/(Ar, Ar). (b) 写出从 x0 到 x1 的计算公式.
(L1 , K1 ) z1 = V1 y1 ∈ K1 r0 − Az1 ⊥ L1 = AK1 (V1T AT AV1 )y1 = V1T AT r0
T (D1 D1 + E T E )y1 = V1T AT r0 T D1 + E T E )−1 V1T AT r0 x1 = x0 + (D1 (1) (1) (2) (1) (1) (2)
A
≤
(κ − 1)2 ek κ2 + 6κ + 1
A,
高等数值分析
2
其中 D1 , D2 是 m×m 的非奇异矩阵. 取 L1 = K1 = Span{e1 , e2 , · · · , em }, L2 = K2 = Span{em+1 , em+2 , · · · , en }. 依次用 (L1 , K1 ), (L2 , K2 ) 按讲义 46 和 47 页 公式 Az ∗ = r0 r0 − Azm ⊥ L W T AV ym = W T r0 xm = x0 + V (W T AV )−1 W T r0 各进行一步计算. 写出一个程序不断按这个方法计算下去, 并验证算法收敛 性. 用 Li = AKi 重复上述各步骤. (0) x1 解. 对任意给定 x0 = (0) , 令 r0 = b − Ax0 , V1 = [e1 , e2 , . . . , em ], V2 = x2 [em , em+1 , . . . , en ]. 对 Li = Ki 情形, 依次用 (L1 , K1 ), (L2 , K2 ) 各进行一步计算： (L1 , K1 ) z1 = V1 y1 r0 − Az1 ⊥ L1 (V1T AV1 )y1 = V1T r0 ,
第 n 步突然收敛。
m−1 r0 } 的 5. 对于 m < n, Arnoldi 过程产生的向量序列 {vi }m i=1 是空间 Span{r0 , Ar0 , · · · , A
一组标准正交基, 而且 vm+1 ⊥ Span{r0 , Ar0 , · · · , Am−1 r0 }.
证明. 第 m(< n) 步的 Arnoldi 过程可表示为
xk = xk−1 + V1 y1 + V2 y2 ; rk = rk−1 − AV1 y1 − AV2 y2 . 收敛性: T (D1 D1 + E T E )−1 rk = r k − 1 − A (I − B )rk−1 算法收敛 ⇔ ρ(I − B ) < 1 ⇔ 0 < λ(B ) < 2.
高等数值分析
0 0 0 1 0 1 0 0 , v1 = . 完成 Arnoldi 过程. 0 1 1 0 1 1 1 0
5
1 0 0 0 0 1/√3 −1/√2 1/√6 解. V = √ √ 0 1/ 3 0 −2/ 6 √ √ √ 0 1/ 3 1/ 2 1/ 6
高等数值分析第三章作业参考答案
1. 考虑线性方程组 Ax = b, 其中 A 是对称正定矩阵. 用 Galerkin 原理求解方 程 K = L = Span(v ), 这里 v 是一个固定的向量. e0 = x∗ − x0 , e1 = x∗ − x1 证明 (e1 , Ae1 ) = (e0 , Ae0 ) − (r, v )2 /(Av, v ), 其中 r = b − Ax0 . v 应当取哪个向量在某种意义上是最佳的? (∗)
(1) (1) (1) (1)
(L2 , K2 ) z1 = V2 y2 r0 − Az1 ⊥ L2 (V2T AV2 )y2 = V2T r0 , −D2 y2 = V2T r0
−1 T x1 = x0 − V2 D2 V2 r0 (2) (2) (2) (2)
D1 y1 = V1T r0
−1 T x1 = x0 + V1 D1 V1 r0
Avm = h1,m v1 + h2,m v2 + · · · + hm,m vm + hm+1,m vm+1 , 其中 v1 = r0 / r0 . 当 m = 1 时, Av1 = h1,1 v1 + h2,1 v2 , 显然命题成立. 设 m = k < n−1 时成立, 即 Span{v1 , v2 , . . . , vk } = Span{r0 , Ar0 , · · · , Ak−1 r0 } 并 且 vk+1 ⊥ Span{r0 , Ar0 , · · · , Ak−1 r0 }. 那么当 m = k + 1 时, 由 Avk = h1,k v1 + h2,k v2 + · · · + hk,k vk + hk+1,k vk+1 , 知 vk+1 ∈ Span{v1 , v2 , . . . , vk , Avk } ⊂ Span{r0 , Ar0 , · · · , Ak r0 }, 因此 Span{v1 , v2 , . . . , vk+1 } ⊂ Span{r0 , Ar0 , · · · , Ak r0 } , 又由于两者都是 k+1 维的, 故 Span{v1 , v2 , . . . , vk+1 } = Span{r0 , Ar0 , · · · , Ak r0 }. 又 vk+2 ⊥ Span{v1 , v2 , . . . , vk+1 } 故 vk+2 ⊥ Span{r0 , Ar0 , · · · , Ak r0 }. 1 1 6. 取 A = 1 1
高等数值分析
3
收敛性: −1 D1 r k = rk −1 − A 0 = −1 ED1 0
−1 −D2
rk −1

−1 −F D2
rk −1
Brk−1
−1 −1 ) < 1. F D2 算法收敛 ⇔ ρ(B ) < 1 ⇔ ρ(ED1
对 Li = AKi 情形, 依次用 (L1 , K1 ), (L2 , K2 ) 各进行一步计算:
rk+1 = rk − βk Ark − γk Apk ; xk+1 = xk + βk rk + γk pk . 此算法本质上是由 CG 迭代一步就重启得到的, 所以是收敛的, 下面给 出证法. 设用此算法得到的 xk+1 = xk + p ¯1 (A)rk , 那么 ek+1
A
=
p1 ∈P1
min ek + p1 (A)rk
证明. 令 x1 = x0 + αv , 那么 r1 = r − αAv, e1 = e0 − αv . 由 Galerkin 原
理, 有 (r1 , v ) = 0, 因 此 α = (r, v )/(Av, v ). 注 意 到 r1 = Ae1 , r = Ae, 有 (Ae1 , v ) = 0. 于是 (e1 , Ae1 ) = (e0 −αv, Ae1 ) = (e0 , Ae1 ) = (e0 , Ae0 )−α(e0 , Av ) = (e0 , Ae0 )−α(r, v ) 即 (∗) 式成立. 由 (∗) 式知当 v = e0 时, e1 到); 在最速下降意义下 v = r 时最佳. 2. 求证： 考虑线性方程组 Ax = b, 其中 A 是对称正定矩阵. 取 K = L = Span(r, Ar). 用 Galerkin 方法求解, 其中 r 是上一步的残余向量. (a) 用 r 和满足 (r, Ap) = 0 的 p 向量构成 K 中的一组基. 给出计算 p 的公 式.
高等数值分析
4

T (D2 D2
+ F F)
T
−1
AT r k − 1
4. 令

0 −1 3 −2 . . . , . . . A= . . . ... ... −2 0 −1 3 xm − x∗
3
−2

1 0 b = . . . 2
解. 设 x1 = x0 + β1 r + β2 p, 那么 r1 = r − β1 Ar − β2 Ap, 再由 Galerkin 原
理, 有 (r1 , r) = (r1 , p) = 0, 解得 β1 = (r, r)/(Ar, r), β2 = (r, p)/(Ap, p).
(c) 该算法收敛吗?
高等数值分析
1
解. 该算法可描述为：
(1) 选初始 x0 ∈ Rn , 计算初始残差 r0 = b − Ax0 , ε > 0 为停机准则； (2) 对 k = 0, 1, 2, . . .直到 rk < ε αk = − (rk , Ark ) ; (Ark , Ark )
pk = rk + αArk ; βk = (rk , rk ) ; (Ark , rk ) (rk , pk ) γk = ; (Apk , pk )
A
A
≤ ek + p ¯1 (A)rk
A
= ek − p ¯1 (A)Aek
A
Biblioteka Baidu
≤ max |p ˜(λi )| ek
1≤i≤n
其中 0 < λ1 ≤ . . . ≤ λn 为 A 的特征值, p ˜(t) = 1 − tp ¯1 (t) 是过 (0, 1) 点的
λn 二次多项式. 当 p ˜ 满足 p ˜(λ1 ) = p ˜(λn ) = −p ˜( λ1 + ) 时可使 max |p ˜(λi )| 达 2 1≤i≤n
到最小. 经计算可得 min max |p ˜(λi )| ≤
p ˜ 1≤i≤n
(λ1 − λn )2 <1 (λ1 − λn )2 + 8λ1 λn
故若令 κ = λ1 /λn , 则 ek+1 方法收敛. 3. 考虑方程组 x b 1 = 1 , −E −D2 x2 b2 D1 −F
1 0 0 0 √ √ 3 2 − 2 / 6 0 H= √ √ 0 2/ 6 1/2 − 3/6 √ 0 0 3/6 −1/2
1 1 7. 取 上 题 中 的 矩 阵 A 以 及 b = , m = 2, x0 = 0, 完 成 求 解 Ax = 1 1 G b 的 Arnoldi 和 GMRES 算法, 得出 xA 2 和 x2 . 8. 分别取 A1 = (aij ) = ( 2 1 ), i+j−1 1
得如下算法： (1) 选初始 x0 ∈ Rn , 计算初始残差 r0 = b − Ax0 , ε > 0 为停机准则； (2) 对 k = 1, 2, . . .直到 rk < ε
T 求解 (D1 D1 + E T E )y1 = (AT rk−1 )(1 : m) ; T 求解 (D2 D2 + F T F )y2 = (AT rk−1 )(m + 1 : n) ;
得如下算法： (1) 选初始 x0 ∈ Rn , 计算初始残差 r0 = b − Ax0 , ε > 0 为停机准则； (2) 对 k = 1, 2, . . .直到 rk < ε 求解 D1 y1 = rk−1 (1 : m) ; 求解 − D2 y2 = rk−1 (m + 1 : n) ; xk = xk−1 + V1 y1 + V2 y2 ; rk = rk−1 − AV1 y1 − AV2 y2 .
用 Galerkin 原理求解 Ax = b. 取 x0 = 0, Vm = Wm = (e1 , e2 , · · · , em ). 对不同的 m, 观察 确解. b − Axm 和 的变化, 其中 x∗ 为方程的精
解. 对于
b − Axm
和
xm − x∗ , 都是前 n − 1 步下降趋势微乎其微, 到
(L2 , K2 ) z1 = V2 y2 ∈ K2 r0 − Az1 ⊥ L2 = AK2 (V2T AT AV2 )y2 = V2T AT r0
T (D2 D2 + F T F )y2 = V2T AT r0 T D2 + F T F )−1 V2T AT r0 x1 = x0 + (D2 (2) (2)