本章主要讲的是求解方程组
AX=b(∗) AX=b\qquad\qquad\qquad\qquad (*) AX=b(∗)
其中 A∈Rn×nA\in R^{n\times n}A∈Rn×n 为非奇异矩阵
Gauss消元法
前提条件
消元过程的所有主元素 akk(k)≠0⇐⇒a_{kk}^{(k)}\neq0\Leftarrow \Rightarrowakk(k)̸=0⇐⇒ 系数矩阵 AAA 的 kkk 阶顺序主子阵 det(Ak)(k=1,2,⋯ ,m)det(A_k)(k=1,2,\cdots,m)det(Ak)(k=1,2,⋯,m) 均非奇异
列选主元
我们从子块(如果是构造上三角矩阵,它的左边全是零)
(ak+1,k+1(k+1)ak+2,k+1(k+1)⋮an,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)⎠⎟⎟⎟⎟⎞
中找到绝对值最大的元素 ap,k+1(k+1)a^{(k+1)}_{p,k+1}ap,k+1(k+1) ,将整个矩阵的第 k+1k+1k+1 行与第 ppp 行互换,从而使每次做消元时,主元素最大。
前推过程
构造形式如下:
A(n)=(a11(1)a12(1)⋯a1n(1)a22(2)⋯a2n(2)⋮ann(n)),b(n)=(b1(1)b1(2)⋮b1(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)⎠⎟⎟⎟⎟⎞
回代过程
我们从第 nnn 个方程开始,自下而上依次解出 xn,xn−1,⋯ ,x1x_n,x_{n-1},\cdots,x_{1}xn,xn−1,⋯,x1 。
Doolittle分解法
我们记
A=LU
A=LU
A=LU
定理: 若矩阵 A∈Rn×nA\in R^{n\times n}A∈Rn×n 的顺序主子式 det(Ai)≠0(i=1,2,⋯ ,n),det(A_i)\neq0(i=1,2,\cdots,n),det(Ai)̸=0(i=1,2,⋯,n), 则存在唯一的下三角矩阵 LLL 及上三角矩阵 UUU 使得上式成立。
求解过程可以分为下列子过程:
LY=b⇒Y=(y1,y2,⋯ ,yn)T⇒UX=Y⇒X.
LY=b\Rightarrow Y=(y_1,y_2,\cdots,y_n)^T\Rightarrow UX=Y\Rightarrow X.
LY=b⇒Y=(y1,y2,⋯,yn)T⇒UX=Y⇒X.
步骤:
- LLL 的第一列与 AAA 的第一列相同;
- 求 UUU 的第一行;
- 求 LLL 的第二列;
- 求 UUU 的第二行;
- ⋯⋯\cdots\cdots⋯⋯
最后可得到 LLL 与 UUU ,在得到解 XXX 。
改进的Cholesky分解法
没看懂,建议直接看《计算方法(第二版)》的P60 。
追赶法
也就是Gauss消元法的特殊应用,没什么难,《计算方法(第二版)》的P62。
扰动分析
条件数 Cond(A):∣∣A−1∣∣∣∣A∣∣Cond(A):||A^{-1}||||A||Cond(A):∣∣A−1∣∣∣∣A∣∣。当 Cond(A)>>1Cond(A)>>1Cond(A)>>1 时,方程组 (∗)(*)(∗) 视为病态的。常用的条件数有:
Cond1(A)=∣∣A−1∣∣1∣∣A∣∣1,Cond∞(A)=∣∣A−1∣∣∞∣∣A∣∣∞.
Cond_1(A)=||A^{-1}||_1||A||_1,\\
Cond_\infty(A)=||A^{-1}||_\infty||A||_\infty.
Cond1(A)=∣∣A−1∣∣1∣∣A∣∣1,Cond∞(A)=∣∣A−1∣∣∞∣∣A∣∣∞.
上述方式就是一般的直接法,而迭代法比直接法更适合于现代大规模科学工程计算。
一般单步迭代法
设线性方程 (∗)(*)(∗) 有如下迭代格式:
X(k+1)=BK(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) 时,迭代格式 (∗∗)(**)(∗∗) 收敛的充要条件是其迭代矩阵 BBB 的谱半径 ρ(B)<1\rho(B)<1ρ(B)<1。
Jacobi迭代法
将线性方程组 (∗)(*)(∗) 的系数矩阵 AAA 分解为
A=L+D+U,
A=L+D+U,
A=L+D+U,
其中 D=diag(a11,a22,⋯ ,ann),D=diag(a_{11},a_{22},\cdots,a_{nn}),D=diag(a11,a22,⋯,ann),
L=(00⋯00a210⋯00a31a32⋯00⋮⋮⋮⋮an1an2⋯an,n−10),
L=\left(\begin{array}{ccccc}
0&0&\cdots&0&0\\
a_{21}&0&\cdots&0&0\\
a_{31}&a_{32}&\cdots&0&0\\
\vdots&\vdots&&\vdots&\vdots\\
a_{n1}&a_{n2}&\cdots&a_{n,n-1}&0\\
\end{array}\right)
,\\
L=⎝⎜⎜⎜⎜⎜⎛0a21a31⋮an100a32⋮an2⋯⋯⋯⋯000⋮an,n−1000⋮0⎠⎟⎟⎟⎟⎟⎞,
U=(0a12a13⋯a1n00a23⋯a2n⋮⋮⋮⋮000⋯an−1,n000⋯0).
U=\left(\begin{array}{ccccc}
0&a_{12}&a_{13}&\cdots&a_{1n}\\
0&0&a_{23}&\cdots&a_{2n}\\
\vdots&\vdots&\vdots&&\vdots\\
0&0&0&\cdots&a_{n-1,n}\\
0&0&0&\cdots&0\\
\end{array}\right)
.\\
U=⎝⎜⎜⎜⎜⎜⎛00⋮00a120⋮00a13a23⋮00⋯⋯⋯⋯a1na2n⋮an−1,n0⎠⎟⎟⎟⎟⎟⎞.
于是有
(L+D+U)X=b⇒DX=−(L+U)X+b⇒X=−D−1(L+D)X+D−1b
(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=b⇒DX=−(L+U)X+b⇒X=−D−1(L+D)X+D−1b
Jacobi迭代公式:
X(k+1)=−D−1(L+D)X(k)+D−1b,k=0,1,⋯ ,
X^{(k+1)}=-D^{-1}(L+D)X^{(k)}+D^{-1}b,\quad k=0,1,\cdots,
X(k+1)=−D−1(L+D)X(k)+D−1b,k=0,1,⋯,
定理(重要): 若线性方程组 (∗)(*)(∗) 的系数矩阵 AAA 严格对角占优,则Jacobi迭代法是收敛的。
Gauss-Seidel迭代法
方程组 (∗)(*)(∗) 也可以等价地写为
(D+L)X=−UX+b
(D+L)X=-UX+b
(D+L)X=−UX+b
类似Jacobi迭代法可以得到Gauss-Seidel迭代法:
X(k+1)=−(D+L)−1UX(k)+(D+L)−1b
X^{(k+1)}=-(D+L)^{-1}UX^{(k)}+(D+L)^{-1}b
X(k+1)=−(D+L)−1UX(k)+(D+L)−1b
定理(重要): 若线性方程组 (∗)(*)(∗) 的系数矩阵 AAA 严格对角占优,则Gauss-Seidel迭代法是收敛的。
JOR迭代法
JOR迭代法是由Jacobi迭代法加入松弛因子 www 得到。
由:
X(k+1)=X(k)+w步长
X^{(k+1)}=X^{(k)}+w步长
X(k+1)=X(k)+w步长
可以得到JOR迭代法:
X(k+1)=X(k)−wD−1(AX(k)−b).
X^{(k+1)}=X^{(k)}-wD^{-1}(AX^{(k)}-b).
X(k+1)=X(k)−wD−1(AX(k)−b).
JOR迭代法有最佳松弛因子
wopt=22−λmaxBJ−λminBJ,
w_{opt}=\frac{2}{2-\lambda^{B_J}_{max}-\lambda^{B_J}_{min}},
wopt=2−λmaxBJ−λminBJ2,
其中 λmaxBJ,λminBJ\lambda^{B_J}_{max},\lambda^{B_J}_{min}λmaxBJ,λminBJ 分别表示Jacobi迭代矩阵 BJ=−D−1(L+U)B_J=-D^{-1}(L+U)BJ=−D−1(L+U) 的最大和最小特征值。此外,当 λmaxBJ≠λminBJ\lambda^{B_J}_{max}\neq\lambda^{B_J}_{min}λmaxBJ̸=λminBJ 时,JOR迭代法的收敛速度相较于对应的Jacobi迭代法的收敛速度快。
定理(重要): 若线性方程组 (∗)(*)(∗) 的系数矩阵 AAA 严格对角占优,则松弛因子 w∈(0,1]w\in (0,1]w∈(0,1] 的JOR迭代法是收敛的。
SOR迭代法
SOR迭代法是由Gauss-Seidel迭代法加入松弛因子 www 得到。
由:
DX(k+1)=DX(k)+w步长
DX^{(k+1)}=DX^{(k)}+w步长
DX(k+1)=DX(k)+w步长
得到SOR迭代法:
X(k+1)=(D+wL)−1{[(1−w)D−wU]X(k)+wb}.
X^{(k+1)}=(D+wL)^{-1}\{[(1-w)D-wU]X^{(k)}+wb\}.
X(k+1)=(D+wL)−1{[(1−w)D−wU]X(k)+wb}.
SOR迭代法的最佳松弛因子
wopt=21+1−ρ2(BJ)
w_{opt}=\frac{2}{1+\sqrt{1-\rho^2(B_J)}}
wopt=1+1−ρ2(BJ)2
定理(重要): 若线性方程组 (∗)(*)(∗) 的系数矩阵 AAA 严格对角占优,则松弛因子 w∈(0,1]w\in (0,1]w∈(0,1] 的SOR迭代法是收敛的。
下面是自己推导Jacobi,Gauss-Seidel,JOR,SOR的过程: