马上注册，结交更多好友，享用更多功能，让你轻松玩转社区。
您需要 登录 才可以下载或查看，没有账号？立即注册
x
系列笔记为上科学计算课上随意记下的笔记，大概有一些缭乱，勉强拼集着看。以记载和梳理为主，证明较少。
一、直接法
1. 高斯消去法
计算复杂度分析：
- 乘除法次数：
∑ k = 1 n [ ( n − k + 2 ) ( n − k ) + ( n − k + 1 ) ] = n 3 3 + n 2 − n 3 \sum_{k=1}^n \left[(n-k+2)(n-k) + (n-k+1)\right] = \frac{n^3}{3} + n^2 - \frac{n}{3} k=1∑n[(n−k+2)(n−k)+(n−k+1)]=3n3+n2−3n
- 加减法次数：
∑ k = 1 n [ ( n − k + 1 ) ( n − k ) ] + ( n − k + 1 ) = n 3 3 + 1 2 n 2 + 1 6 n \sum_{k=1}^n \left[(n-k+1)(n-k)\right] +(n-k+1)= \frac{n^3}{3} +\frac{1}{2}n^{2}+\frac{1}{6}n k=1∑n[(n−k+1)(n−k)]+(n−k+1)=3n3+21n2+61n
- 总复杂度： O ( n 3 ) O(n^3) O(n3)，以 n 3 3 \frac{n^3}{3} 3n3 为主项
2. 选主元消去法
- 列主元法：在列中选择最大元素作为主元，避免小主元导致的数值不稳定
- 全主元法：在整个子矩阵中选择最大主元（计算量更大，但稳定性最佳）
3. LU 分解
定理：若矩阵 A ∈ R n × n A \in \mathbb{R}^{n \times n} A∈Rn×n 的次序主子式均非零，则存在唯一的单元下三角矩阵 L L L 和上三角矩阵 U U U，使得 A = L U A = LU A=LU。
存储优化：
- U U U 存储于原矩阵的上三角部分
- L L L 的非单元元素存储于下三角部分（对角元1无需存储）
Crout 分解：将 A A A 分解为下三角矩阵与单元上三角矩阵的乘积。
明白LU分解并推得递推表达式，可以从初等矩阵乘积进行消元的角度进行明白，也可以直接从LU的形式入手，直接进行计算。从初等矩阵乘积消元的角度看，利用列主元LU分解，也就是选取列主元消元，我们可以得到以下结论：
列主元 LU 分解
定理：对恣意非奇异矩阵 A A A，存在置换矩阵 P P P，使得：
P A = L U PA = LU PA=LU
此中 L L L 为单元下三角矩阵， U U U 为上三角矩阵。
实现步骤：
假设前r－1行都进行了换行消元的过程，即 l i k l_{ik} lik ( i = 1 , ⋯ , n ; k = (i=1,\cdots,n;k= (i=1,⋯,n;k= 1 , ⋯ r − 1 ) 1,\cdots r-1) 1,⋯r−1) ）及 u k i u_{ki} uki ( i = 1 , ⋯ , n ; k = 1 , ⋯ r − 1 ) (\boldsymbol{i}=1,\cdots,n;k=1,\cdots r-1) (i=1,⋯,n;k=1,⋯r−1) 均已得到。下面要得到必要互换的两行以及U的第r行。换行的尺度是希望形成的 u r r u_{rr} urr 绝对值最大。考虑A的第r列的r至n元素 a i r a_{ir} air ,i=r…n,有
a i r = ∑ k = 1 r − 1 l i k u k r + u r r ^ a_{ir}=\sum_{k=1}^{r-1}l_{ik}u_{kr}+\widehat{u_{rr}} air=k=1∑r−1likukr+urr
此中 s i = s_{i}= si= u r r ^ \widehat{u_{rr}} urr 为各行若称为主行后，大概得到的 u r r u_{rr} urr ．比较 ∣ s i r ∣ = max ∣ s i ∣ \left|s_{i_r}\right|=\max|s_i| ∣sir∣=max∣si∣ ，互换第r行与 i r i_r ir 行，
4. 对称矩阵的特殊分解
1. LDLᵀ 分解
定理：若对称矩阵 A A A 的次序主子式非零，则存在唯一分解：
A = L D L T A = LDL^T A=LDLT
此中 L L L 为单元下三角阵， D D D 为对角阵。
2. Cholesky 分解
定理：若 A A A 对称正定，则存在唯一的下三角矩阵 L L L（对角元全正），使得：
A = L L T A = LL^T A=LLT
计算步骤：
l k k = a k k − ∑ m = 1 k − 1 l k m 2 l i k = 1 l k k ( a i k − ∑ m = 1 k − 1 l i m l k m ) ( i > k ) \begin{aligned} l_{kk} &= \sqrt{a_{kk} - \sum_{m=1}^{k-1} l_{km}^2} \\ l_{ik} &= \frac{1}{l_{kk}} \left( a_{ik} - \sum_{m=1}^{k-1} l_{im}l_{km} \right) \quad (i > k) \end{aligned} lkklik=akk−m=1∑k−1lkm2 =lkk1(aik−m=1∑k−1limlkm)(i>k)
二、迭代法
1. 根本理论
定理1：下列命题等价
- lim k → ∞ B k = 0 \lim_{k\to\infty}B^k=0 limk→∞Bk=0
- ρ ( B ) < 1 \rho(B)<1 ρ(B)<1
- 存在 ∥ ⋅ ∥ ε \|{\cdot}\|_{\varepsilon} ∥⋅∥ε，使 ∥ B ∥ ε < 1 \left\|\boldsymbol{B}\right\|_{\varepsilon}<\mathbf{1} ∥B∥ε<1
定理2：对恣意一种范数，有
lim k → ∞ ∣ ∣ B k ∣ ∣ 1 / k = ρ ( B ) \lim_{k\to\infty}\left|\left|B^{k}\right|\right|^{1/k}=\rho(B) k→∞lim Bk 1/k=ρ(B)
这个定理告诉我们，无论我们选取任何一种范数，其收敛举动都是由谱半径决定的。
下面先容迭代法的根本定理。
迭代格式： x ( k + 1 ) = B x ( k ) + g x^{(k+1)} = Bx^{(k)} + g x(k+1)=Bx(k)+g
收敛定理：迭代法收敛当且仅当 ρ ( B ) < 1 \rho(B) < 1 ρ(B)<1（谱半径小于1）
毛病估计：
∥ x ∗ − x ( k ) ∥ ≤ ∥ B ∥ k ⋅ ∥ x ( ∗ ) − x ( 0 ) ∥ ∥ x ∗ − x ( k ) ∥ ≤ ∥ B ∥ k 1 − ∥ B ∥ ∥ x ( 1 ) − x ( 0 ) ∥ ∥ x ∗ − x ( k ) ∥ ≤ ∥ B ∥ 1 − ∥ B ∥ ∥ x ( k ) − x ( k − 1 ) ∥ \begin{aligned} \|x^* - x^{(k)}\| &\leq \|B\|^k\cdot \|x^{(*)} - x^{(0)}\| \\ \|x^* - x^{(k)}\| &\leq \frac{\|B\|^k}{1 - \|B\|} \|x^{(1)} - x^{(0)}\| \\ \|x^* - x^{(k)}\| &\leq \frac{\|B\|}{1 - \|B\|} \|x^{(k)} - x^{(k-1)}\| \end{aligned} ∥x∗−x(k)∥∥x∗−x(k)∥∥x∗−x(k)∥≤∥B∥k⋅∥x(∗)−x(0)∥≤1−∥B∥∥B∥k∥x(1)−x(0)∥≤1−∥B∥∥B∥∥x(k)−x(k−1)∥
收敛速度分析
对恣意给定的范数，对恣意 δ > 0 \delta>0 δ>0 ,存在着正整数K使得 k ≥ K k\geq K k≥K 时，
∣ ∣ ε ( k ) ∣ ∣ = ∣ ∣ B k ε ( 0 ) ∣ ∣ ≤ ∣ ∣ B k ∣ ∣ ∣ ∣ ε ( 0 ) ∣ ∣ ≤ ( ρ ( B ) + δ ) k ∣ ∣ ε ( 0 ) ∣ ∣ \left|\left|\varepsilon^{(k)}\right|\right|=\left|\left|B^{k}\varepsilon^{(0)}\right|\right|\leq\left|\left|B^{k}\right|\right|\left|\left|\varepsilon^{(0)}\right|\right|\leq\left(\rho(B)+\delta\right)^{k}\left|\left|\varepsilon^{(0)}\right|\right| ε(k) = Bkε(0) ≤ Bk ε(0) ≤(ρ(B)+δ)k ε(0)
于是有:
∣ ∣ ε ( k ) ∣ ∣ ∣ ∣ ε ( 0 ) ∣ ∣ ≤ ( ρ ( B ) + δ ) k ≤ τ \frac{||\varepsilon^{(k)}||}{||\varepsilon^{(0)}||}\leq(\rho(B)+\delta)^k\leq\tau ∣∣ε(0)∣∣∣∣ε(k)∣∣≤(ρ(B)+δ)k≤τ
于是为达到精度 τ \tau τ，需迭代步数：
k ≥ − ln τ − ln ( ρ ( B ) + δ ) k \geq \frac{-\ln \tau}{-\ln (\rho(B) + \delta)} k≥−ln(ρ(B)+δ)−lnτ
2. 矩阵分裂法
将 A A A 分裂为 A = M − N A = M - N A=M−N，迭代公式：
x ( k + 1 ) = M − 1 N x ( k ) + M − 1 b x^{(k+1)} = M^{-1}Nx^{(k)} + M^{-1}b x(k+1)=M−1Nx(k)+M−1b
Richardson迭代
将A分裂为：
A = I − ( I − A ) A = \mathcal{I} - (\mathcal{I} - A) A=I−(I−A)
于是：
x ( k ) = ( I − A ) x ( k − 1 ) + b x^{(k)} = (\mathcal{I}-A)x^{(k-1)}+b x(k)=(I−A)x(k−1)+b
可以写成更有启发性的形式：
x ( k ) = x ( k − 1 ) + ( b − A x ( k − 1 ) ) x^{(k)} = x^{(k-1)} + \left({b-Ax^{(k-1)}}\right) x(k)=x(k−1)+(b−Ax(k−1))
以及更一般的形式：
x ( k ) = x ( k − 1 ) + α ( b − A x ( k − 1 ) ) x^{(k)} = x^{(k-1)} + \alpha\left({b-Ax^{(k-1)}}\right) x(k)=x(k−1)+α(b−Ax(k−1))
一般来说，在Richardson迭代中，我们要求矩阵 A A A是一个对称正定矩阵，此中正定矩阵是因为要求谱半径绝对值要求严酷小于1。于是迭代矩阵的特征值：
λ [ I − α A ] = 1 − α A \lambda[\mathcal{I}-\alpha A] = 1-\alpha A λ[I−αA]=1−αA
收敛条件为：
α λ max ( A ) < 2 \alpha \lambda_{\max\: }(A) <2 αλmax(A)<2
外推法：
x ( k ) = ν [ ( I − A ) x ( k − 1 ) + b ] + ( 1 − ν ) x ( k − 1 ) 充要条件: 0 < ν < 2 λ m a x . \begin{aligned}&x^{(k)}\:=\:\nu\Big[(I-A)x^{(k-1)}\:+\:b\Big]+(1-\nu)x^{(k-1)}\\&\text{充要条件:}0<\nu<\frac{2}{\lambda_{max}}.\end{aligned} x(k)=ν[(I−A)x(k−1)+b]+(1−ν)x(k−1)充要条件:0<ν<λmax2.
Jacobi迭代与Gauss-Seidel迭代
- Jacobi 迭代： M = D M = D M=D（对角阵）
x i ( k ) = 1 a i i ( b i − ∑ j ≠ i a i j x j ( k − 1 ) ) x_i^{(k)} = \frac{1}{a_{ii}} \left( b_i - \sum_{j \neq i} a_{ij}x_j^{(k-1)} \right) xi(k)=aii1 bi−j=i∑aijxj(k−1)
- 可并行计算，需存储两代迭代值
- 存储需求：需保存两个向量 x ( k ) , x ( k + 1 ) x^{(k)},x^{(k+1)} x(k),x(k+1)。
- 收敛速度：通常较慢，但对角占优或对称正定矩阵包管收敛。
- Gauss-Seidel 迭代： M = D − L M = D - L M=D−L（下三角部分）
x i ( k ) = 1 a i i ( b i − ∑ j < i a i j x j ( k ) − ∑ j > i a i j x j ( k − 1 ) ) x_i^{(k)} = \frac{1}{a_{ii}} \left( b_i - \sum_{j < i} a_{ij}x_j^{(k)} - \sum_{j > i} a_{ij}x_j^{(k-1)} \right) xi(k)=aii1(bi−j<i∑aijxj(k)−j>i∑aijxj(k−1))
收敛性分析
设 A = ( a i j ) n × n A=\left(a_{ij}\right)_{n\times n} A=(aij)n×n
如果满足
∣ a i i ∣ > ∑ j = 1 , j ≠ i n a i j ( i = 1 , 2 , ⋯ , n ) \begin{aligned}|a_{ii}|&>\sum_{j=1,j\neq i}^{n}a_{ij}\quad(i=1,2,\cdots,n)\end{aligned} ∣aii∣>j=1,j=i∑naij(i=1,2,⋯,n)
即每一行的对角线元素绝对值都严酷大于同行其他元素的绝对值之和， 则称A为严酷对角占优阵。如果满足
∣ a i i ∣ ≥ ∑ j = 1 j ≠ i n a i j ( i = 1 , 2 , ⋯ , n ) |a_{ii}|\geq\sum_{j=1\atop j\neq i}^{n}a_{ij}\quad(i=1,2,\cdots,n) ∣aii∣≥j=ij=1∑naij(i=1,2,⋯,n)
且至少有一个不等式严酷建立，则称A为弱对角占优阵。如果存在 n n n 阶置换矩阵 P P P 使
P T A P = [ A 11 A 12 0 A 22 ] P^TAP=\begin{bmatrix}A_{11}&A_{12}\\0&A_{22}\end{bmatrix} PTAP=[A110A12A22]
建立，此中 A 11 A_{11} A11 为r阶子矩阵， A 22 A_{22} A22 为 n − r n-r n−r 阶子矩阵，则称A 是可约矩阵。否则，称A为不可约矩阵。
定理1：如果 A = ( a i j ) A=(a_{ij} ) A=(aij)是严酷对角占优矩阵或不可约弱对角占优矩阵，则Jacobi和Gauss-Seidel方法均收敛
起首分析Jacobi方法。等价于证明矩阵：
B J = I − D − 1 A B_{J} = \mathcal{I}-D^{-1}A BJ=I−D−1A
的所有特征值的模都小于1。等价于证明对于恣意模大于1的特征值 λ \lambda λ：
∣ λ I − B J ∣ ≠ 0 \lvert \lambda \mathcal{I} - B_J\rvert\neq 0 ∣λI−BJ∣=0
而：
∣ λ I − B J ∣ = ∣ λ I − D − 1 ( L + U ) ∣ = ∣ D − 1 ∣ ⋅ ∣ λ D − ( L + U ) ∣ \begin{align*} \lvert \lambda \mathcal{I}-B_{J}\rvert &= \lvert \lambda \mathcal{I} - D^{-1}\left({L+U}\right)\rvert\\ & = \lvert D^{-1}\rvert\cdot\lvert \lambda D - \left({L+U}\right)\rvert \end{align*} ∣λI−BJ∣=∣λI−D−1(L+U)∣=∣D−1∣⋅∣λD−(L+U)∣
于是这是显然的。
定理2：如果 A A A 对称， a i i > 0 a_{ii}>0 aii>0则：
1.Gauss-Seidel方法收敛的充分须要条件是A正定；
2.Jacobi方法收敛的充分须要条件是A与2D-A均正定。
这里我们只证明Jacobi方法。
B J = I − D − 1 A ⟹ λ ( B J ) = λ ( I − D − 1 A ) = 1 − λ ( D − 1 A ) B_{J} = \mathcal{I} - D^{-1}A\implies\lambda \left({B_{J}}\right) = \lambda \left({\mathcal{I} - D^{-1}A}\right) = 1 - \lambda \left({D ^{-1}A}\right) BJ=I−D−1A⟹λ(BJ)=λ(I−D−1A)=1−λ(D−1A)
而：
λ ( D − 1 A ) = λ ( D − 1 / 2 A D − 1 / 2 ) \lambda(D^{-1}A) = \lambda \left({D^{-1/2}AD^{-1/2}}\right) λ(D−1A)=λ(D−1/2AD−1/2)
于是：
λ ( D − 1 A ) > 0 ⟺ A 正定 \lambda(D^{-1}A)>0 \iff A正定 λ(D−1A)>0⟺A正定
λ ( D − 1 A ) = λ ( D − 1 / 2 A D − 1 / 2 ) < 2 ⟺ 2 D − A 正定 \lambda(D^{-1}A) = \lambda(D^{-1/2}AD^{-1/2})<2 \iff 2\mathcal{D} - A正定 λ(D−1A)=λ(D−1/2AD−1/2)<2⟺2D−A正定
SOR 迭代
在Gauss-Seidel迭代的底子上引入松懈因子 ω \omega ω，通过调整 ω \omega ω加速收敛。 迭代格式为Gauss-Seidel效果的加权平均：
x i ( k ) = ( 1 − ω ) x i ( k − 1 ) + ω a i i ( b i − ∑ j = 1 i − 1 a i j x j ( k ) − ∑ j = i + 1 n a i j x j ( k − 1 ) ) x_i^{(k)} = (1-\omega)x_i^{(k-1)} + \frac{\omega}{a_{ii}} \left( b_i - \sum_{j=1}^{i-1} a_{ij}x_j^{(k)} - \sum_{j=i+1}^n a_{ij}x_j^{(k-1)} \right) xi(k)=(1−ω)xi(k−1)+aiiω(bi−j=1∑i−1aijxj(k)−j=i+1∑naijxj(k−1))
其矩阵分裂形式如下
- 分裂矩阵：
A = 1 ω D − L − ( 1 − ω ω D + U ) A = \frac{1}{\omega}D - L - \left( \frac{1-\omega}{\omega}D + U \right) A=ω1D−L−(ω1−ωD+U)
于是对应迭代矩阵：
B SOR = ( D − ω L ) − 1 [ ( 1 − ω ) D + ω U ] B_{\text{SOR}} = (D - \omega L)^{-1} \left[ (1-\omega)D + \omega U \right] BSOR=(D−ωL)−1[(1−ω)D+ωU]
收敛性定理
- 须要条件：
0 < ω < 2 0 < \omega < 2 0<ω<2
证明可以通过计算:
det B S O R = 1 det D ⋅ det D ( 1 − w ) n = ( 1 − w ) n \det B_{SOR} = \frac{1}{\det D}\cdot\det D\left({1-w}\right)^{n} = \left({1-w}\right)^{n} detBSOR=detD1⋅detD(1−w)n=(1−w)n
如果 w ∉ ( 0 , 2 ) w \not \in (0,2) w∈(0,2)，则行列式必然大于1，于是必然存在大于1的特征值，从而不收敛。但是也应该明白如果行列式小于1，也不能包管收敛。
- 对称正定矩阵的收敛性（Ostrowski定理）：
若 A A A对称正定且 0 < ω < 2 0 < \omega < 2 0<ω<2，则SOR方法收敛。
这个定理的证明比较复杂，这里略去。
- 最优松懈因子：
当 A A A为对称正定的三对角矩阵时，最优松懈因子为：
ω opt = 2 1 + 1 − ρ ( B J ) 2 \omega_{\text{opt}} = \frac{2}{1 + \sqrt{1 - \rho(B_J)^2}} ωopt=1+1−ρ(BJ)2 2
此中 ρ ( B J ) \rho(B_J) ρ(BJ)为Jacobi迭代矩阵的谱半径。
免责声明：如果侵犯了您的权益，请联系站长，我们会及时删除侵权内容，谢谢合作！更多信息从访问主页：qidao123.com:ToB企服之家，中国第一个企服评测及商务社交产业平台。 |
