学习笔记十六——线性方程组的数值解法

本章主要讲的是求解方程组
A X = b ( ∗ ) AX=b\qquad\qquad\qquad\qquad (*) AX=b()
其中 A ∈ R n × n A\in R^{n\times n} ARn×n 为非奇异矩阵

Gauss消元法

前提条件

消元过程的所有主元素 a k k ( k ) ≠ 0 ⇐ ⇒ a_{kk}^{(k)}\neq0\Leftarrow \Rightarrow akk(k)̸=0 系数矩阵 A A A k k k 阶顺序主子阵 d e t ( A k ) ( k = 1 , 2 , ⋯   , m ) det(A_k)(k=1,2,\cdots,m) det(Ak)(k=1,2,,m) 均非奇异

列选主元

我们从子块(如果是构造上三角矩阵,它的左边全是零)
( a k + 1 , k + 1 ( k + 1 ) a k + 2 , k + 1 ( k + 1 ) ⋮ a n , k + 1 ( k + 1 ) ) \left(\begin{array}{ccccc} a^{(k+1)}_{k+1,k+1}\\ a^{(k+1)}_{k+2,k+1}\\ \vdots\\ a^{(k+1)}_{n,k+1} \end{array}\right) ak+1,k+1(k+1)ak+2,k+1(k+1)an,k+1(k+1)
中找到绝对值最大的元素 a p , k + 1 ( k + 1 ) a^{(k+1)}_{p,k+1} ap,k+1(k+1) ,将整个矩阵的第 k + 1 k+1 k+1 行与第 p p p 行互换,从而使每次做消元时,主元素最大。

前推过程

构造形式如下:
A ( n ) = ( a 11 ( 1 ) a 12 ( 1 ) ⋯ a 1 n ( 1 ) a 22 ( 2 ) ⋯ a 2 n ( 2 ) ⋮ a n n ( n ) ) , b ( n ) = ( b 1 ( 1 ) b 1 ( 2 ) ⋮ b 1 ( n ) ) A^{(n)}=\left(\begin{array}{ccccc} a^{(1)}_{11}&a^{(1)}_{12}&\cdots&a^{(1)}_{1n}\\ &a^{(2)}_{22}&\cdots&a^{(2)}_{2n}\\ &&&\vdots\\ &&&a^{(n)}_{nn}\\ \end{array}\right) ,\quad b^{(n)}= \left(\begin{array}{ccccc} b_1^{(1)}\\ b_1^{(2)}\\ \vdots\\ b_1^{(n)} \end{array}\right) A(n)=a11(1)a12(1)a22(2)a1n(1)a2n(2)ann(n),b(n)=b1(1)b1(2)b1(n)

回代过程

我们从第 n n n 个方程开始,自下而上依次解出 x n , x n − 1 , ⋯   , x 1 x_n,x_{n-1},\cdots,x_{1} xn,xn1,,x1

Doolittle分解法

我们记
A = L U A=LU A=LU
定理: 若矩阵 A ∈ R n × n A\in R^{n\times n} ARn×n 的顺序主子式 d e t ( A i ) ≠ 0 ( i = 1 , 2 , ⋯   , n ) , det(A_i)\neq0(i=1,2,\cdots,n), det(Ai)̸=0(i=1,2,,n), 则存在唯一的下三角矩阵 L L L 及上三角矩阵 U U U 使得上式成立。
求解过程可以分为下列子过程:
L Y = b ⇒ Y = ( y 1 , y 2 , ⋯   , y n ) T ⇒ U X = Y ⇒ X . LY=b\Rightarrow Y=(y_1,y_2,\cdots,y_n)^T\Rightarrow UX=Y\Rightarrow X. LY=bY=(y1,y2,,yn)TUX=YX.
步骤:

  1. L L L 的第一列与 A A A 的第一列相同;
  2. U U U 的第一行;
  3. L L L 的第二列;
  4. U U U 的第二行;
  5. ⋯ ⋯ \cdots\cdots

最后可得到 L L L U U U ,在得到解 X X X

改进的Cholesky分解法

没看懂,建议直接看《计算方法(第二版)》的P60 。

追赶法

也就是Gauss消元法的特殊应用,没什么难,《计算方法(第二版)》的P62。

扰动分析

条件数 C o n d ( A ) : ∣ ∣ A − 1 ∣ ∣ ∣ ∣ A ∣ ∣ Cond(A):||A^{-1}||||A|| Cond(A):A1A。当 C o n d ( A ) > > 1 Cond(A)>>1 Cond(A)>>1 时,方程组 ( ∗ ) (*) () 视为病态的。常用的条件数有:
C o n d 1 ( A ) = ∣ ∣ A − 1 ∣ ∣ 1 ∣ ∣ A ∣ ∣ 1 , C o n d ∞ ( A ) = ∣ ∣ A − 1 ∣ ∣ ∞ ∣ ∣ A ∣ ∣ ∞ . Cond_1(A)=||A^{-1}||_1||A||_1,\\ Cond_\infty(A)=||A^{-1}||_\infty||A||_\infty. Cond1(A)=A11A1,Cond(A)=A1A.
上述方式就是一般的直接法,而迭代法比直接法更适合于现代大规模科学工程计算。

一般单步迭代法

设线性方程 ( ∗ ) (*) () 有如下迭代格式:
X ( k + 1 ) = B K ( k ) + F , k = 0 , 1 , 2 , ⋯   , ( ∗ ∗ ) X^{(k+1)}=BK^{(k)}+F,\quad k=0,1,2,\cdots,\qquad(**) X(k+1)=BK(k)+F,k=0,1,2,,()
定理(重要): 当给定初始向量 X ( 0 ) X^{(0)} X(0) 时,迭代格式 ( ∗ ∗ ) (**) () 收敛的充要条件是其迭代矩阵 B B B 的谱半径 ρ ( B ) &lt; 1 \rho(B)&lt;1 ρ(B)<1

Jacobi迭代法

将线性方程组 ( ∗ ) (*) () 的系数矩阵 A A A 分解为
A = L + D + U , A=L+D+U, A=L+D+U,
其中 D = d i a g ( a 11 , a 22 , ⋯ &ThinSpace; , a n n ) , D=diag(a_{11},a_{22},\cdots,a_{nn}), D=diag(a11,a22,,ann),
L = ( 0 0 ⋯ 0 0 a 21 0 ⋯ 0 0 a 31 a 32 ⋯ 0 0 ⋮ ⋮ ⋮ ⋮ a n 1 a n 2 ⋯ a n , n − 1 0 ) , L=\left(\begin{array}{ccccc} 0&amp;0&amp;\cdots&amp;0&amp;0\\ a_{21}&amp;0&amp;\cdots&amp;0&amp;0\\ a_{31}&amp;a_{32}&amp;\cdots&amp;0&amp;0\\ \vdots&amp;\vdots&amp;&amp;\vdots&amp;\vdots\\ a_{n1}&amp;a_{n2}&amp;\cdots&amp;a_{n,n-1}&amp;0\\ \end{array}\right) ,\\ L=0a21a31an100a32an2000an,n10000,
U = ( 0 a 12 a 13 ⋯ a 1 n 0 0 a 23 ⋯ a 2 n ⋮ ⋮ ⋮ ⋮ 0 0 0 ⋯ a n − 1 , n 0 0 0 ⋯ 0 ) . U=\left(\begin{array}{ccccc} 0&amp;a_{12}&amp;a_{13}&amp;\cdots&amp;a_{1n}\\ 0&amp;0&amp;a_{23}&amp;\cdots&amp;a_{2n}\\ \vdots&amp;\vdots&amp;\vdots&amp;&amp;\vdots\\ 0&amp;0&amp;0&amp;\cdots&amp;a_{n-1,n}\\ 0&amp;0&amp;0&amp;\cdots&amp;0\\ \end{array}\right) .\\ U=0000a12000a13a2300a1na2nan1,n0.
于是有
( L + D + U ) X = b ⇒ D X = − ( L + U ) X + b ⇒ X = − D − 1 ( L + D ) X + D − 1 b (L+D+U)X=b\\ \Rightarrow DX=-(L+U)X+b\\ \Rightarrow X=-D^{-1}(L+D)X+D^{-1}b (L+D+U)X=bDX=(L+U)X+bX=D1(L+D)X+D1b
Jacobi迭代公式:
X ( k + 1 ) = − D − 1 ( L + D ) X ( k ) + D − 1 b , k = 0 , 1 , ⋯ &ThinSpace; , X^{(k+1)}=-D^{-1}(L+D)X^{(k)}+D^{-1}b,\quad k=0,1,\cdots, X(k+1)=D1(L+D)X(k)+D1b,k=0,1,,
定理(重要): 若线性方程组 ( ∗ ) (*) () 的系数矩阵 A A A 严格对角占优,则Jacobi迭代法是收敛的。

Gauss-Seidel迭代法

方程组 ( ∗ ) (*) () 也可以等价地写为
( D + L ) X = − U X + b (D+L)X=-UX+b (D+L)X=UX+b
类似Jacobi迭代法可以得到Gauss-Seidel迭代法:
X ( k + 1 ) = − ( D + L ) − 1 U X ( k ) + ( D + L ) − 1 b X^{(k+1)}=-(D+L)^{-1}UX^{(k)}+(D+L)^{-1}b X(k+1)=(D+L)1UX(k)+(D+L)1b
定理(重要): 若线性方程组 ( ∗ ) (*) () 的系数矩阵 A A A 严格对角占优,则Gauss-Seidel迭代法是收敛的。

JOR迭代法

JOR迭代法是由Jacobi迭代法加入松弛因子 w w w 得到。
由:
X ( k + 1 ) = X ( k ) + w 步 长 X^{(k+1)}=X^{(k)}+w步长 X(k+1)=X(k)+w
可以得到JOR迭代法:
X ( k + 1 ) = X ( k ) − w D − 1 ( A X ( k ) − b ) . X^{(k+1)}=X^{(k)}-wD^{-1}(AX^{(k)}-b). X(k+1)=X(k)wD1(AX(k)b).
JOR迭代法有最佳松弛因子
w o p t = 2 2 − λ m a x B J − λ m i n B J , w_{opt}=\frac{2}{2-\lambda^{B_J}_{max}-\lambda^{B_J}_{min}}, wopt=2λmaxBJλminBJ2,
其中 λ m a x B J , λ m i n B J \lambda^{B_J}_{max},\lambda^{B_J}_{min} λmaxBJ,λminBJ 分别表示Jacobi迭代矩阵 B J = − D − 1 ( L + U ) B_J=-D^{-1}(L+U) BJ=D1(L+U) 的最大和最小特征值。此外,当 λ m a x B J ≠ λ m i n B J \lambda^{B_J}_{max}\neq\lambda^{B_J}_{min} λmaxBJ̸=λminBJ 时,JOR迭代法的收敛速度相较于对应的Jacobi迭代法的收敛速度快。
定理(重要): 若线性方程组 ( ∗ ) (*) () 的系数矩阵 A A A 严格对角占优,则松弛因子 w ∈ ( 0 , 1 ] w\in (0,1] w(0,1] 的JOR迭代法是收敛的。

SOR迭代法

SOR迭代法是由Gauss-Seidel迭代法加入松弛因子 w w w 得到。
由:
D X ( k + 1 ) = D X ( k ) + w 步 长 DX^{(k+1)}=DX^{(k)}+w步长 DX(k+1)=DX(k)+w
得到SOR迭代法:
X ( k + 1 ) = ( D + w L ) − 1 { [ ( 1 − w ) D − w U ] X ( k ) + w b } . X^{(k+1)}=(D+wL)^{-1}\{[(1-w)D-wU]X^{(k)}+wb\}. X(k+1)=(D+wL)1{[(1w)DwU]X(k)+wb}.
SOR迭代法的最佳松弛因子
w o p t = 2 1 + 1 − ρ 2 ( B J ) w_{opt}=\frac{2}{1+\sqrt{1-\rho^2(B_J)}} wopt=1+1ρ2(BJ) 2
定理(重要): 若线性方程组 ( ∗ ) (*) () 的系数矩阵 A A A 严格对角占优,则松弛因子 w ∈ ( 0 , 1 ] w\in (0,1] w(0,1] 的SOR迭代法是收敛的。

下面是自己推导Jacobi,Gauss-Seidel,JOR,SOR的过程:
在这里插入图片描述

  • 0
    点赞
  • 2
    收藏
    觉得还不错? 一键收藏
  • 0
    评论

“相关推荐”对你有帮助么?

  • 非常没帮助
  • 没帮助
  • 一般
  • 有帮助
  • 非常有帮助
提交
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值