数值分析复习笔记-第六章-线性方程组的迭代解法

Chapter6 线性方程组的迭代解法

6.1 范数

6.1.1 向量范数

1. 定义：对任意n维向量$\mathbf{x}$，有对应的一个非负实数$||\mathbf{x}||$满足->
   • 非负性：$||x||⩾0$，并且$||x||=0\to x=0$
   • 齐次性：对于任意实数$\gamma$，有$||\gamma x||=\gamma ||x||$
   • 三角不等性：$||x||+||y||⩾||x+y||$
2. 分类
   • 1-norms: norm(x,1)
     $$||v||=|v_1|+|v_2|+|v_3|+\cdots +|v_n|=\sum _{i=1}^n|v_i|$$
   • 2-norms（欧几里得范数）: norm(x,2)
     $$\Vert x{\Vert }_2=\sqrt{v_1^2+v_2^2+\cdots +v_n^2}=\sqrt{\sum _{i=1}^n{v}_i^2}$$
   • $\infty$-norms: norm(x,inf)
     $$\Vert x{\Vert }_{\infty }=\underset{1⩽j⩽n}{\max }\left|x_j\right|$$

6.1.2 矩阵范数

1. 定义：对任意n阶方阵$\mathbf{A}$，有对应的一个非负实数$||\mathbf{A}||$满足->
   • 非负性：$||\mathbf{A}||⩾0$，并且$||\mathbf{A}||=0\to \mathbf{A}=0$
   • 齐次性：对于任意实数$\gamma$，有$||\gamma \mathbf{A}||=\gamma ||\mathbf{A}||$
   • 三角不等性：$||\mathbf{A}||+||\mathbf{B}||⩾||\mathbf{A}+\mathbf{B}||$
   • 相容性：$||\mathbf{AB}||\le ||\mathbf{A}||\cdot ||\mathbf{B}||$
2. 谱半径：矩阵特征值的最大值
   $$\rho (\mathbf{A})=\underset{1\le j\le n}{\max }|{\lambda }_i|$$
3. 分类
   • 1-norms(列模): norm(A,1)
     $$\Vert A{\Vert }_1=\underset{1\le j\le n}{\max }\sum _{i=1}^n\left|a_{ij}\right|$$
   • 2-norms(谱模): norm(A,2)
     $$\Vert A{\Vert }_2=\sqrt{\rho (A^{T}A)}$$
   • $\infty$-norms(行模): norm(A,inf)
     $$\Vert A{\Vert }_1=\underset{1\le i\le n}{\max }\sum _{j=1}^n\left|a_{ij}\right|$$
   • F-norms(vector 2-norms)： norm(A,'fro')
     $$\Vert \mathrm{A}{\Vert }_{\mathrm{F}}={\left(\sum _{\mathrm{i}=1}^{\mathrm{m}}\sum _{\mathrm{j}=1}^{\mathrm{n}}{\left|{\mathrm{a}}_{\mathrm{ij}}\right|}^2\right)}^{\frac{1}{2}}$$

6.1.3 相容关系

1. 定义：若如下不等式成立，我们称向量范数和矩阵范数相容
   $$||\mathbf{Ax}||\le ||\mathbf{A}||\cdot ||\mathbf{x}||$$
2. 其余范数应满足：
   $$\begin{gathered} ||\mathbf{Ax}|{|}_1\le ||\mathbf{A}|{|}_1\cdot ||\mathbf{x}|{|}_1 \\ ||\mathbf{Ax}|{|}_2\le ||\mathbf{A}|{|}_2\cdot ||\mathbf{x}|{|}_2 \\ ||\mathbf{Ax}|{|}_{\infty }\le ||\mathbf{A}|{|}_{\infty }\cdot ||\mathbf{x}|{|}_{\infty } \\ ||\mathbf{Ax}|{|}_2\le ||\mathbf{A}|{|}_F\cdot ||\mathbf{x}|{|}_2 \end{gathered}$$

6.2 条件数

6.2.1 条件数的推导

 对一般的非奇异线性方程组$\mathbf{Ax}=\mathbf{b}$，在$\mathbf{x}$和$\mathbf{b}$两端分别增加扰动$\mathbf{\delta b}$和$\mathbf{\delta x}$，有：
$$A(x+\delta x)=b+\delta b$$
 带入$\mathbf{Ax}=\mathbf{b}$，得到：
$$\delta x=A^{-1}\delta b$$
 两边同时取范数,根据相容性得：
$$\Vert \delta x\Vert \le \Vert A^{-1}\Vert \cdot \Vert \delta b\Vert$$
 对$\mathbf{Ax}=\mathbf{b}$两边取范数，得到：
$$\Vert b\Vert \le \Vert A\Vert \cdot \Vert x\Vert$$
$$\frac{1}{\Vert x\Vert }\le \Vert A\Vert \cdot \frac{1}{\Vert b\Vert }$$
 两边分别相乘：
$$\frac{\Vert \delta x\Vert }{\Vert x\Vert }\le \Vert A\Vert \Vert A^{-1}\Vert \frac{\Vert \delta b\Vert }{\Vert b\Vert }$$
式中这个$\Vert A\Vert \Vert A^{-1}\Vert$玩意起到一个放大倍数的作用，对方程组得解的相对误差起到关键性作用。

6.2.2 定义

若$\mathbf{A}$为n阶非奇异矩阵，称$cond(\mathbf{A})=||\mathbf{A}||\cdot ||{\mathbf{A}}^{-1}||$为条件数

6.2.3 分类

常用的有：

1. cond(A,inf)：
   $$con{d}_{\infty }(A)=||A|{|}_{\infty }\cdot ||A^{-1}|{|}_{\infty }$$
2. cond(A,2)：
   $$con{d}_2(A)=||A|{|}_2\cdot ||A^{-1}|{|}_2$$
3. 如果A是对称正定矩阵，有:
   $$con{d}_2(A)=\frac{{\lambda }_1(A)}{{\lambda }_n(A)}=\frac{A特征值最大值}{A特征值最小值}$$
   证明：
   $$Ax=\lambda x=>AAx=A\lambda x=\lambda Ax={\lambda }^{2}x$$
   $$A\cdot A^T=A^2=>if{\lambda }_1^2>...>{\lambda }_n^2$$
   $$\begin{gathered} ||A\cdot A^T|{|}_2=\Vert A{\Vert }_2=\sqrt{\rho (A^{T}A)}={\lambda }_1 \\ ||A^{-1}\cdot (A^T{)}^{-1}|{|}_2=\Vert A{\Vert }_2=\sqrt{\rho (A^{T}A)}=\frac{1}{{\lambda }_n} \end{gathered}$$

6.2.4 性质

1. 对于任意的n阶非奇异矩阵A，$cond(A)\ge 1$
   $$||A||\cdot ||A^{-1}||\ge ||A\cdot A^{-1}||=||I||=1$$
2. 对于任意的n阶非奇异矩阵A和常数c，$cond(cA)=cond(A)$
3. 对于任意正交矩阵P，$con{d}_2(P)=1$，且$con{d}_2(AP)=con{d}_2(PA)=con{d}_2(A)$

   选用常数和正交矩阵，无法改善条件数

6.2.5 作用

 条件数越小，表示这个矩阵式良态矩阵，反之，为病态矩阵。

6.3 基本迭代法

定常迭代法的迭代矩阵通常保持不变

6.3.1 basic concept

1. 本质：find a matrix, derive the origin matrix to a new form
   $$\begin{gathered} origin:Ax=b \\ \downarrow \\ A=M-n \\ Mx=Nx+b,M=非奇异 \\ x=Bx+g=\{\begin{array}{ll}B& =M^{-1}N\\ g& =M^{-1}b\end{array} \\ \downarrow \\ newform:x^{(k+1)}=B\cdot x^{(k)}+g \end{gathered}$$
2. 定义：我们希望$M$和$N$具有特殊性质，因此，定义对角矩阵$D$、上三角矩阵$L$和下三角矩阵$U$：
   $$\begin{gathered} A=D-L-U \\ D=diag(a_{11},a_{22},...,a_{nn}) \\ L=-\left(\begin{array}{cccc}0& & & \\ a_{21}& 0& & \\ ⋮& & \ddots & \\ a_{n1}& a_{n2}& \cdots & 0\end{array}\right) \\ U=-\left(\begin{array}{cccc}0& a_{12}& \cdots & a_{1n}\\ & 0& \cdots & a_{2n}\\ & & \ddots & ⋮\\ & & & 0\end{array}\right) \end{gathered}$$

6.3.2 Jacobi

1. 本质：
   $$x^{(k+1)}=D^{-1}(L+U)x^{(k)}+D^{-1}b$$
2. 推导：
   $$\{\begin{array}{l}{a}_{11}{x}_1+a_{12}{x}_2+\cdots a_{1n}{x}_n=b_1\\ a_{21}{x}_1+a_{22}{x}_2+\cdots a_{2n}{x}_n=b_2\\ \cdots \\ a_{n1}{x}_1+a_{n2}{x}_2+\cdots a_{nn}{x}_n=b_n\end{array}\iff x=Mx+g$$
   $$\begin{gathered} \Downarrow \\ \{\begin{array}{l}{x}_1=-\frac{a_{12}}{a_{11}}{x}_2-\cdots -\frac{a_{1n}}{a_{11}}{x}_n+\frac{b_1}{a_{11}}\\ x_2=-\frac{a_{21}}{a_{22}}{x}_1-\cdots -\frac{a_{2n}}{a_{22}}{x}_n+\frac{b_2}{a_{22}}\\ \cdots \\ x_n=-\frac{a_{n1}}{a_{nn}}{x}_1-\frac{a_{n2}}{a_{nn}}{x}_2-\cdots -\frac{a_{nn-1}}{a_{nn}}{x}_{n-1}+\frac{b_n}{a_{nn}}\end{array} \end{gathered}$$
   $$\begin{gathered} \Downarrow \\ {x}_i^{(k+1)}=\frac{b_i-{\sum }_{j=1}^{i-1}{a}_{ij}{x}_j^{(k)}-{\sum }_{j=i+1}^n{a}_{ij}{x}_j^{(k)}}{a_{ii}};i=1,2,\cdots ,n \end{gathered}$$
3. 例题一道

写一个Jacobi迭代程序，输入维数n，求解Ax=b，其中：
$$U=-\left(\begin{array}{cccc}n+1& 1& \cdots & 1\\ 1& n+2& \cdots & 1\\ ⋮& & \ddots & ⋮\\ 1& 1& \cdots & 2n\end{array}\right),b=\left(\begin{array}{c}1\\ 2\\ ⋮\\ n\end{array}\right)$$

function [x,iteration]=jacobi(n,tol)
% n为维数，tol为误差
A=ones(n,n);
D=diag([n+1:2*n]); %n*n
L=-tril(A,-1); %n*n
U=-triu(A,1); % n*n
A=D-L-U;
b=[1:n]'; %n*1
x=zeros(size(b)) % initialize x
for iteration=1:2000
    x=D\(L*x+U*x+b) % remember divide by'\' no '/'
    error=norm(b-A*x,2)/norm(b,2) % 2-norm to calculate error
    if error<tol
        break;
    end
end
end

4. notes
   简单记忆：向左倒是左除\，向右倒是右除/
   /：右除：a/b表示矩阵a乘以矩阵b的逆
   \：左除：a\b表示矩阵a的逆乘以b。

6.3.3 Gauss-Seidel

1. 本质：
   $$x^{(k+1)}=(D-L{)}^{-1}U{x}^{(k)}+(D-L{)}^{-1}b$$
2. 上一次迭代出来的x，会继承给下一次迭代
   $$a_{ii}\cdot x_i^{(k+1)}+\sum _{j=1}^{i-1}{a}_{ij}{x}_j^{(k+1)}=\left(b_i-\sum _{j=i+1}^n{a}_{ij}{x}_j^{(k)}\right),i=1,2,\dots ,n$$
   $$\begin{gathered} \Downarrow \\ (D-L)\cdot x=Ux+b \end{gathered}$$
   $$\begin{gathered} \Downarrow \\ {x}_i^{(k+1)}=\frac{1}{a_{ii}}\left(b_i-\sum _{j=1}^{i-1}{a}_{ij}{x}_j^{(k+1)}-\sum _{j=i+1}^n{a}_{ij}{x}_j^{(k)}\right),i=1,2,\dots ,n \end{gathered}$$
3. 一道例题

用高斯赛德尔迭代求解Ax=b，其中：
$$A=\left(\begin{array}{ccc}20& 4& 6\\ 4& 20& 8\\ 6& 8& 20\end{array}\right),b=\left(\begin{array}{c}10\\ -24\\ -22\end{array}\right)$$

A=[20 4 6;4 20 8;6 8 20];
b=[10 -24 -22]';
[xgs,itergx]=gs(A,b,5e-5)
%% guess-seidel
function [x,iter]=gs(A,b,tol)
D=diag(diag(A));
L=-tril(A,-1);
U=-triu(A,1);
x=zeros(size(b));
for iter=1:2000
    x=(D-L)\(U*x+b); % D->(D-L)
    error=norm(b-A*x)/norm(b);
    if error<tol
        break;
    end
end
end

6.3.4 SOR：Successive Over Relaxation

1. 本质：上一次迭代结果+残差修正量
   $$x^{(k+1)}=x^{(k)}+\omega \cdot D^{-1}(L{x}^{(k+1)}+U{x}^{(k)}-D{x}^{(k)}+b)$$
2. 推导：高斯赛德尔迭代->分离$x^{(K)}$，然后对剩余项乘以一个松弛因子
   $$x_i^{(k+1)}=\frac{1}{a_{ii}}\left(b_i-\sum _{j=1}^{i-1}{a}_{ij}{x}_j^{(k+1)}-\sum _{j=i+1}^n{a}_{ij}{x}_j^{(k)}\right),i=1,2,\dots ,n$$
   $$\begin{gathered} \Downarrow \\ {\mathrm{x}}_{\mathrm{i}}^{(\mathrm{k}+1)}={\mathrm{x}}_{\mathrm{i}}^{(\mathrm{k})}+\frac{1}{{\mathrm{a}}_{\mathrm{ii}}}{\mathrm{r}}_{\mathrm{i}}^{(\mathrm{k})} \\ {\mathrm{r}}_{\mathrm{i}}^{(\mathrm{k})}=\left({\mathrm{b}}_{\mathrm{i}}-\sum _{\mathrm{j}=1}^{\mathrm{i}-1}{\mathrm{a}}_{\mathrm{ij}}{\mathrm{x}}_{\mathrm{j}}^{(\mathrm{k}+1)}-\sum _{\mathrm{j}=\mathrm{i}}^{\mathrm{n}}{\mathrm{a}}_{\mathrm{ij}}{\mathrm{x}}_{\mathrm{j}}^{(\mathrm{k})}\right),\mathrm{i}=1,2,\cdots ,\mathrm{n} \end{gathered}$$
   $$\begin{gathered} \Downarrow \\ {\mathrm{x}}_{\mathrm{i}}^{(\mathrm{k}+1)}=(1-\omega ){\mathrm{x}}_{\mathrm{i}}^{(\mathrm{k})}+\frac{\omega }{{\mathrm{a}}_{\mathrm{ii}}}\left({\mathrm{b}}_{\mathrm{i}}-\sum _{\mathrm{j}=1}^{\mathrm{i}-1}{\mathrm{a}}_{\mathrm{ij}}{\mathrm{x}}_{\mathrm{j}}^{(\mathrm{k}+1)}-\sum _{\mathrm{j}=\mathrm{i}+1}^{\mathrm{n}}{\mathrm{a}}_{\mathrm{ij}}{\mathrm{x}}_{\mathrm{j}}^{(\mathrm{k})}\right),\mathrm{i}=1,2,\cdots ,\mathrm{n}(2.5) \end{gathered}$$
3. 一道例题

用超松弛迭代求解Ax=b，其中：
$$A=\left(\begin{array}{ccc}20& 4& 6\\ 4& 20& 8\\ 6& 8& 20\end{array}\right),b=\left(\begin{array}{c}10\\ -24\\ -22\end{array}\right)$$

function [x,iter]=sor(A,b,omega,tol)
D=diag(diag(A));
L=-tril(A,-1);
U=-triu(A,1);
x=zeros(size(b));
for iter=1:2000
    x=(D-omega*L)\((1-omega)*D*x+omega*U*x+omega*b);
    error=norm(b-A*x)/norm(b);
    if error<tol
        break;
    end
end
end

4. 特点
   • 高斯赛德尔迭代是$\omega =1$时的SOR迭代
   • 采用松弛因子是为了加快迭代方法收敛的速度

6.4 迭代法收敛性和误差分析

6.4.1 可约矩阵

1. 排列矩阵：
   • 每行每列仅有唯一的非0元1的方阵P
   • 对任意矩阵A左乘一个排列矩阵P，就是对矩阵A的行进行排列。
   • 对任意矩阵A右乘一个排列矩阵P，就是对矩阵A的列进行排列。
     $$P=\left(\begin{array}{ccc}0& 1& 0\\ 1& 0& 0\\ 0& 0& 1\end{array}\right)$$
2. 可约矩阵：设A是n阶矩阵，如果存在n阶排列矩阵P，满足：
   $$P^{T}AP=\left(\begin{array}{cc}{A}_{11}& A_{12}\\ 0& A_{22}\end{array}\right)$$
   其中，$A_{11}$和$A_{22}$分别为r阶和n-r阶的方阵，则称为可约矩阵；若不存在，则称为不可约矩阵。
3. 应用：类似于高斯消元，将x和b也进行类似转换，可以得到，从而加快运算
   $$\{\begin{array}{l}{A}_{11}{y}_1+A_{12}{y}_2=f_1\\ A_{11}{y}_2=f_2\end{array}$$

6.4.2 对角占优

1. 若n阶矩阵A满足：
   $$\left|a_{ii}\right|>\sum _{j\ne i}\left|a_{ij}\right|,i=1,\cdots ,n.$$
   • 至少满足一个：弱对角占优
   • 全部满足：严格对角占优
   • 都不满足：不对角占优
2. 例题一道

判断下列矩阵属于什么占优?
$$\left(\begin{array}{ccc}5& 3& 2\\ -2& 4& 2\\ 6& 1& 8\end{array}\right)$$
5=|3|+|2|;
4=|-2|+|2|;
8>|6|+|1|;
因此，这是一个弱对角占优

3. 引理一：若A是严格对角占优OR不可约弱对角占优，则A是可逆矩阵（非奇异矩阵==det|A|≠0）

反证法，若$det|A|=0$，则对$Ax=0$有非零解$x=(x_1,\cdots ,x_n)$
令$x_r$为最大的一个，则$|x_r|>0$，不妨令$|x_r|=1$，取第r个方程：
$$\sum _{j=1}^n{a}_{rj}{x}_j=0$$
$$-a_{rr}{x}_r=\sum _{j=1,j\ne r}^n{a}_{rj}{x}_j$$
$$|a_{rr}||x_r|\le \sum _{j=1,j\ne r}^n|a_{rj}||x_j|$$
$$|a_{rr}|\le \sum _{j=1,j\ne r}^n|a_{rj}||x_j|\le \sum _{j=1,j\ne r}^n|a_{rj}|$$
因此，A可逆

4. 引理二：n阶矩阵A的k次幂${\mathbf{A}}^{\mathbf{k}}\to 0(\mathbf{k}\to \infty )$ $\iff$ 谱半径$\mathbf{\rho }(\mathbf{A})<1$

6.4.3 迭代法基本定理

1. 定理一：$\mathbf{\rho }(\mathbf{A})\le ||\mathbf{A}||$，其中，A表示任一与某一向量范数相容的矩阵范数
2. 定理二：对于基本迭代格式，给定初值$x^{(0)}$，$x^{\ast }$是真值，有下列收敛结果和误差估计
   • 迭代格式收敛的充要条件：$\mathbf{\rho }(\mathbf{B})<1$
   • 如果$||B||<1$，则有：
     $$\begin{array}{rl}& \Vert x^{(k)}-x_{\ast }\Vert ⩽\frac{\Vert B{\Vert }^k}{1-\Vert B\Vert }\Vert x^{(1)}-x^{(0)}\Vert \\ & \Vert x^{(k)}-x_{\ast }\Vert ⩽\frac{\Vert B\Vert }{1-\Vert B\Vert }\Vert {}^{(k)}-x^{(k-1)}\Vert \end{array}$$

proof1:
$$\begin{gathered} x^{(k)}-x^{\ast }=B(x^{(k-1)}-x^{\ast })=\cdots =B^k(x^{(0)}-x^{\ast }) \\ \downarrow \\ {x}^{(0)}-x^{\ast }=const \\ \downarrow \\ 由引理二：whenk\to \infty ,B^k\to 0\leftrightarrow \rho (B)<1 \end{gathered}$$
proof2: 运用矩阵范数的三角不等式原理
$$\begin{gathered} \begin{array}{rl}{x}^{(k)}-x^{\ast }& =x^{(k)}-x^{(k+1)}+x^{(k+1)}-x^{\ast }\\ & =[(B{x}^{(k-1)}+f)-(B{x}^{(k)}+f)]+[(B{x}^{(k)}+f)-(B{x}^{\ast }+f)]\\ & =B\left(x^{(k-1)}-x^{(k)}\right)+B\left(x^{(k)}-x^{\ast }\right)\end{array} \\ \downarrow \\ \begin{array}{rl}\Vert x^{(k)}-x^{(k+1)}+x^{(k+1)}-x^{\ast }\Vert & ⩽\Vert B{x}^{(k-1)}+f-\left(B{x}^{(k)}+f\right)\Vert +\Vert \left(B{x}^{(k)}+f\right)-\left(B{x}^{\ast }+f\right)\Vert \\ & ⩽\Vert B\Vert \Vert x^{(k-1)}-x^{(k)}\Vert +\Vert B\Vert \Vert x^{(k)}-x^{\ast }\Vert \end{array} \\ \downarrow \\ [1-\Vert B\Vert )\Vert x^{(k)}-x^{\ast }\Vert ⩽\Vert B\Vert \Vert x^{(k-1)}-x^{(k)}\Vert \\ \downarrow \\ \Vert x^{(k)}-x^A\Vert ⩽\frac{\Vert B\Vert }{(1-\Vert B\Vert )}\Vert x^{(k-1)}-x^{(k)}\Vert \end{gathered}$$

1. 定理三：若A是严格对角占优OR不可约弱对角占优 $\Rightarrow$ Jacobi和GS迭代法收敛
2. 定理四：若A是对称正定矩阵，Jacobi收敛 $\iff$ 2D-A也为正定对称
3. 定理五：SOR迭代收敛 $\Rightarrow$ $0<\mathbf{\omega }<2$

$B_{SOR}=(D-\omega L{)}^{-1}[(1-\omega )D+\omega U]$，其特征值为${\lambda }_1,\cdots ,{\lambda }_n$
n个谱半径肯定大于各特征值之积，各特征值的积就是det（A）
$$\begin{array}{rl}(\rho (B_{SOR}))& \ge |{\lambda }_1{\lambda }_2\cdots {\lambda }_n|=det(B_{SOR})\\ & =det|下三角矩阵|\cdot det|上三角矩阵|\\ & =|a_{11}^{-1}\cdots a_{nn}^{-1}||(1-\omega {)}^n{a}_{11}\cdots a_{nn}|\\ & =(1-\omega {)}^n<1\end{array}$$

1. 定理六：若A是对称正定矩阵，且$0<\mathbf{\omega }<2$ $\iff$ SOR迭代收敛
2. 一道例题

设线性方程组Ax=b的系数矩阵为：
$$A=\left(\begin{array}{ccc}1& a& a\\ a& 1& a\\ a& a& 1\end{array}\right)$$
证明：
(1) 当-0.5<a<1时，用GS迭代求解收敛
(2) 当-0.5<a<0.5时，用Jacobi迭代求解收敛
(3) 当a=0.8时，用Jacobi迭代求解发散
求解：
(1)由于不是对角占优，因此采用定理六，判断对称正定矩阵，$\omega =1$
(2)定理三，充分条件，证明严格对角占优矩阵
(3)定理四，充要条件，证明2D-A不是对称正定矩阵

6.5 不定常迭代

6.5.1 特点

基于变分方法来最小化线性方程组的残量方法，没有明显的迭代矩阵

6.5.2 分类

• 最速下降法：求解对称正定线性方程组
• 共轭梯度法：求解对称正定线性方程组，本质上是一种变分的方法，即求二次函数的极值
• 广义极小残量法：求解不对称线性方程组

e4468069-f8a3-4ca5-afbf-462bd10c1997