数值分析读书笔记（2）求解线性代数方程组的直接方法

数值分析读书笔记（2）求解线性代数方程组的直接方法

1.引言

矩阵的数值计算一般可以分为直接法和间接法

本章主要介绍Ax=b这类线性方程组求解的直接法，数值求解该方程组的基础思想是Gauss消元法

实质是通过一组满秩的初等行变换，将A保秩变换成一个三角矩阵U，此变换过程称为矩阵A的非奇异上三角化

我们的目的就是寻求一个矩阵P，使得PA=U，其中U是一个三角矩阵，其中Ax=b和Ux=\overline b同解（\overline b=Pb）,有效的生成一个P是我们主要研究的问题

2.初等下三角矩阵--Guass变换矩阵

回顾一下线性代数中的三个初等线性变换

  • 数乘
  • 倍加
  • 互换

我们引入一个一般意义上的初等变换矩阵，它把许多常用的线性变换统一在一个框架里面，在数值线性代数中起着重要的意义

Def： 称C^{n \times n}中如下形式的矩阵E(u,v;\sigma)为初等矩阵:
E(u,v;\sigma)=I-\sigma u v^{H}
其中非零向量u,v\in C^{n},\sigma \neq 0是实或者复数，即
E(u,v;\sigma)=\begin{pmatrix} 1-\sigma u_1v_1 & -\sigma u_1v_2& \cdots &-\sigma u_1v_n\\ -\sigma u_2v_1& 1-\sigma u_2v_2 & \cdots & -\sigma u_2v_n\\ \vdots & \vdots &\ddots &\vdots\\ -\sigma u_nv_1&-\sigma u_nv_2 &\cdots &1-\sigma u_nv_n \end{pmatrix}
选取不同的u，v，\sigma，可以得到许多常用的线性变换矩阵

  • 数乘（E_1=E(e_i,e_i;1-\alpha)）
  • 倍加（E_2=E(e_i,e_j;-\mu)）
  • 互换（E_3=E(e_i-e_j,e_i-e_j;1)）

下面引出初等变换矩阵的一些重要的数学性质
1.两相同向量u，v组成的初等变换矩阵可交换，其积仍然为一个初等矩阵
E(u,v;\sigma)E(u,v;\tau)=E(u,v;\sigma+\tau-\sigma \tau v^{H}u)

证明：
\begin{align} E(u,v;\sigma)E(u,v;\tau) &=(I-\sigma uv^{H}))(I-\tau uv^{H})\\ &=I-\sigma uv^{H}-\tau uv^{H}+ \end{align}

2.若1-\sigma v^{H} u \neq 0 ,则初等矩阵E(u,v;\sigma)可逆，其逆矩阵也是初等矩阵
E^{-1}(u,v;\sigma)=E(u,v;\tau),\tau=\frac{\sigma}{\sigma v^{H}u-1}

3.设v^{\bot}表示和v正交的（n-1）维子空间
a.若u\notin v^{\bot},则E(u,v;\sigma)有n个线性无关的特征向量，该组特征向量由u和v^{\bot}中任取一组基向量组成
b.若u\in v^{\bot},则E(u,v;\sigma)仅有n-1个线性无关的特征向量，该组特征向量由v^{\bot}中任取一组基向量组成

4.det(E(u,v;\sigma))=1-\sigma v^{H}u

5.对任意非零向量a,b\in C^{n},必可适当选取u,v,\sigma使得E(u,v;\sigma)a=b

事实上只需要取u,v,\sigma满足v^{H}a\neq 0,\sigma u=\frac{a-b}{v^{H}a}

由初等变换矩阵引出Guass变换矩阵，我们选取\sigma =-1,u=l_{k}=(0,\cdots,0,l_{k+1},\cdots,l_{nk})^{T},v=e_{k}=(0,\cdots,0,1,0,\cdots,0)^{T},k=1,2,\cdots,n-1

得到n-1个Guass变换矩阵
L_{k}(l_{k})\equiv E(l_{k},e_{k};-1)=I+l_{k}e_{k}=\begin{pmatrix} 1\\& \ddots\\&&1\\&&l_{k+1,k}&1\\&&l_{k+2,k}&&1\\&&\vdots&&&\ddots\\&&l_

下面给出Guass变换矩阵的一些性质

1.det(L_{k})=1+e_{k}^{T}l_{k}=1

2.Guass变换矩阵的逆只需要将\sigma从-1变成+1

3.L_{1}(l_{1})L_{2}(l_{2})\cdots L_{n-1}(l_{n-1})=\begin{pmatrix}1\\l_{21}&1\\l_{31}&l_{32}&1\\\vdots&\vdots&\vdots&\ddots\\l_{n1}&l_{n2}&l_{n3}&\cdots&1\end{pmatrix}

注意左乘的顺序

3.Gauss消元法

先介绍一下顺序Gauss消元法，大概分两步

  • 消元过程
  • 回代过程

在消元过程中，我们不断去左乘Gauss变换矩阵，不断将原矩阵的下三角部分一列列变成0，从而最终变换成一个上三角矩阵

需要注意的是，在一列列的消元过程中，我们需保证a_{ii}\neq 0 (i=1,2,\dots,n),所以需要利用行互换来保证此条件

当然这一切消元过程的前提是，矩阵A应该是非奇异的

经过n-1次的Gauss消元，我们可以得到一个上三角矩阵
L_{n-1}^{-1}\cdots L_{k}^{-1}\cdots L_{2}^{-1}L_{1}^{-1}A^{(1)}=A^{(n)}x=L_{n-1}^{-1}\cdots L_{k}^{-1}\cdots L_{2}^{-1}L_{1}^{-1}b^{(1)}=b^{(n)}

在回代过程中，由于我们得到了一个上三角矩阵，那么就可以从最底行开始逐步解出x

Gauss消元法的复杂度是O(n^{3})，高阶状态下比起克拉默法则运算量要小得多

Gauss消元法过程中，在对各列进行消元的时候，如果主元比较小的话，运算的结果会产生较大的误差，故引入Gauss列主元消元法，即在每一次利用主元消元的步骤之前，把该列中绝对值最大的数所在的行与主元所在的行进行交换

4.三角分解法

我们利用Gauss变换矩阵对Gauss消元法进行进一步的分析
L_{n-1}^{-1}\cdots L_{k}^{-1}\cdots L_{2}^{-1}L_{1}^{-1}A^{(1)}x=A^{(n)}x=Ux
故\begin{align}A&=L_{1}\cdots L_{k}\cdots L_{n-2}L_{n-1}U\\&=(I+l_{1}e_{1}^{T})\cdots(I+l_{n-2}e_{n-2}^{T})(I+l_{n-1}e_{n-1}^{T})U\\&=(I+l_{1}e_{1}^{T}+\cdots+l_{n-1}e_{n-1}^{T})U\end{align}

由此引出矩阵的LU分解，又称Doolittle分解
A=\begin{pmatrix}1\\l_{21}&1\\l_{31}&l_{32}&1\\\vdots&\vdots&\vdots&\ddots\\l_{n1}&l_{n2}&l_{n3}&\cdots&1\end{pmatrix}\begin{pmatrix}u_{11}&u_{12}&u_{13}&\cdots &u_{1n}\\&u_{22}&u_{23}&\cdots&u_{2n}\\&&u_{33}&\cdots &u_{3n}\\&&&\ddots &\vdots \\&&&&u_{nn}\end{pmatrix}

这里再介绍一下Crout分解，即A=LU中的L是一个下三角矩阵，U是单位上三角矩阵

注意到某些特殊矩阵的三角分解也是比较特殊的，这里引入一类带状对角形矩阵
A=\begin{pmatrix}a_{11}&\cdots &a_{1,s+1}&\\\vdots&\ddots&&\ddots\\a_{r+1,1}&&a_{r+1,s+1}&&a_{n-s,n}\\&\ddots&&\ddots\\&&a_{n,n-r}&&a_{nn}\end{pmatrix}

上半带宽为s，下半带宽为r，存在LU分解，其中L是下半带宽为r的单位下三角矩阵，U是上半带宽为s的上三角矩阵

对于r=s=1的这一类更加特殊的矩阵，称为三对角矩阵，对于此类矩阵的三角分解，介绍一种“追赶法”

首先做Crout分解
A=LU=\begin{pmatrix}p_{1} &&&&0\\r_{2}&p_{2}\\&\ddots&\ddots\\&&r_{n-1}&p_{n-1}\\0&&&r_{n}&p_{n}\end{pmatrix}\begin{pmatrix}1&q_{1}&&&0\\&1&q_{2}\\&&\ddots&\ddots\\&&&1&q_{n-1}\\0&&&&1\end{pmatrix}
然后分两步解决此类问题
追：解Ly=b
赶：解Ux=y

注意到正定对称矩阵的三角分解也是特殊的，这里引入Cholesky分解

首先利用Doolittle分解，得A=LU,对U进一步提取对角矩阵diag(u_{11},\dots,u_{nn}),从而有
U=DD^{-1}U=D(D^{-1}U)=DU_{0}
故，A=LDU_{0}，由于A对称正定，A^{T}=A,所以有
A^{T}=(LDU_{0})^{T}=U_{0}^{T}DL^{T}=A=LDU_{0}
由于分解的唯一性，可知L^{T}=U_{0}，从而有
A=LDL^{T}
我们可以记，D^{1/2}=diag(\sqrt{u_{11}},\dots,\sqrt{u_{nn}}),从而
A=LD^{1/2}D^{1/2}L^{T}=LD^{1/2}(LD^{1/2})^{T}=L_{1}L_{1}^{T}
此种分解手段称为Cholesky分解，限定对角元素为正，此类分解唯一

上述的Cholesky分解中涉及了开方的运算，下面介绍一种改进的平方根法

易知，A=LDL^{T},则Ax=LDL^{T}x

先解Ly=b,后解L^{T}x=D^{-1}y,其中D的逆只需要将对角元素取倒数即可

5.向量和矩阵的范数

范数是比长度更为一般的概念，有了范数就可以更好的去测度误差的大小

关于向量范数

Def: V是数域R/C上的线性空间，对于V中任意的元素x，\\if 存在一个唯一的实函数N（x）与之对应，记为\begin{Vmatrix}x\end{Vmatrix},而且需满足三个条件\\1.非负正定，2.齐次性，3.三角不等式

对于非负正定，当仅当x=0，有N（x）=0，否则N（x）> 0;

对于齐次性，有\begin{Vmatrix}\alpha x\end{Vmatrix}=\begin{vmatrix}\alpha\end{vmatrix} \begin{Vmatrix}x\end{Vmatrix} ,\alpha \in K

对于三角不等式，有\begin{Vmatrix}x+y\end{Vmatrix}\leq\begin{Vmatrix}x\end{Vmatrix}+\begin{Vmatrix}y\end{Vmatrix},\forall x,y \in V

这里介绍几种常见的向量范数

  • l_{1}-范数 向量中的元素的绝对值之和
  • l_{2}-范数 向量中的元素的绝对值的平方加起来然后开方
  • l_{\infty}-范数 向量元素中的最大绝对值（使用Cauchy-Schwarz不等式证明三角不等式）
  • l_{p}-范数 向量中的元素的绝对值的p次方加起来然后开p次方根（利用赫尔德不等式即可证明三角不等式）

在最优化理论中可能会涉及加权范数，A为对称正定矩阵，(x^{T}Ax)^{1/2}是一种向量范数，记为\begin{Vmatrix}x\end{Vmatrix}_{A}

在无限维线性空间中，比如在[a,b]区间中，对于所有的实连续函数集合C[a,b],对于其中的一个元素f(x)也是有类似定义的范数

  • 1范数
    \begin{Vmatrix}f(x)\end{Vmatrix}_1=\int _{a}^{b}\begin{vmatrix}f(x)\end{vmatrix}dx
  • p范数
    \begin{Vmatrix}f(x)\end{Vmatrix}_p=\left(\int _{a}^{b}\begin{vmatrix}f(x)\end{vmatrix}^{p}dx \right)^{\frac{1}{p}}
  • ∞范数
    \begin{Vmatrix}f(x)\end{Vmatrix}_\infty=max\begin{vmatrix}f(x)\end{vmatrix} ,a\leq x\leq b

下面介绍一下范数的等价性

对于任意两个定义好的范数，存在两个与向量x无关的非零正常数c1，c2，有
c_1\begin{Vmatrix}x\end{Vmatrix}_\alpha \leq \begin{Vmatrix}x\end{Vmatrix}_\beta\leq c_2\begin{Vmatrix}x\end{Vmatrix}_\alpha
称两个范数等价

不难验证，此处的等价性满足数学定义中的等价性的三个条件，即自反，对称，传递

关于矩阵范数

矩阵范数不仅仅满足非负正定，齐次和三角不等式，而且须满足矩阵相乘的相容性，即
\begin{Vmatrix}AB\end{Vmatrix} \leq \begin{Vmatrix}A\end{Vmatrix}\begin{Vmatrix}B\end{Vmatrix}

这里给出一类特殊的范数， Frobenius范数
\begin{Vmatrix}A\end{Vmatrix}_F=\left(\sum_{j=1}^{m}\sum_{i=0}^{n}\begin{vmatrix}a_{ij}\end{vmatrix}^{2}\right)^{\frac{1}{2}}

对于C^{m \times n}上面的任意一种向量诱导范数，都有\begin{Vmatrix} I\end{Vmatrix} = \max \limits_{\begin{Vmatrix} x\end{Vmatrix}=1} \{ \begin{Vmatrix} Ix\end{Vmatrix}=1 \}

这里给出一种范数的定义，即诱导矩阵范数，诱导矩阵范数和向量范数密切相关

定义：设在两个向量空间C^{m},C^{n}中存在向量范数\begin{Vmatrix} \bullet\end{Vmatrix}_{V}, 定义在C^{m \times n}空间上的矩阵A的由向量范数 \begin{Vmatrix}\bullet\end{Vmatrix}_{V}诱导所给出的矩阵范数为（其中x不为零向量）\begin{Vmatrix} A\end{Vmatrix}_{V}=max\frac{\begin{Vmatrix} Ax\end{Vmatrix}_{V}}{\begin{Vmatrix} x\end{Vmatrix}_{V}}

我们为了解决这个最大值的问题，继续等价定义来优化这个问题

\begin{Vmatrix} A\end{Vmatrix}_{V}=max\frac{\begin{Vmatrix} Ax\end{Vmatrix}_{V}}{\begin{Vmatrix} x\end{Vmatrix}_{V}}=max\begin{Vmatrix} Ax\end{Vmatrix}_{V}
其中第一个max条件为x不为零向量，第二个max条件为\begin{Vmatrix} x\end{Vmatrix}_{V}=1

我们利用诱导范数的定义可以从原来的向量范数中诱导出三种范数，分别是

1范数：对矩阵的每一列中的元素取绝对值之后求和，然后选取其中的最大列作为1范数
2范数：矩阵的最大奇异值，也就是矩阵与矩阵的转置的乘积的最大特征值
无穷范数：对于矩阵的每一行的元素取绝对值之后求和，然后选取其中的最大行作为无穷范数

关于矩阵的应用，这里引入一个Banach引理

设矩阵A属于n*m的复矩阵空间，对于该空间上的某种矩阵范数\begin{Vmatrix} \bullet \end{Vmatrix}_{V},有\begin{Vmatrix} A\end{Vmatrix}_{V}<1,则矩阵（I±A）非奇异，且有\begin{Vmatrix} (I-A)^{-1}\end{Vmatrix}_{V}\leq \frac{\begin{Vmatrix} I\end{Vmatrix}}{1-\begin{Vmatrix} A\end{Vmatrix}}

给出矩阵谱半径的定义

矩阵的谱半径为矩阵的最大特征值，关于矩阵的谱半径，它不超过其任意一种矩阵范数（当矩阵是Hermite矩阵时，矩阵的2范数恰好等于矩阵的谱半径）

继续给出线性方程组中条件数的定义

在某一矩阵空间中，对于某一矩阵范数，矩阵的条件数=矩阵的范数×矩阵的逆的范数，即
Cond(A)=\begin{Vmatrix} A\end{Vmatrix}_{V} \times \begin{Vmatrix} A^{-1}\end{Vmatrix}_{V}

对于矩阵的条件数来说，它显然大于等于1，当矩阵恰好是正交矩阵的时候，矩阵的条件数恰好等于1
当矩阵为对称阵，对应的矩阵范数为2范数的时候，此时的条件数称之为谱条件数，其值等于最大特征值除以最小特征值，然后取绝对值

©著作权归作者所有,转载或内容合作请联系作者
  • 人面猴
    序言：七十年代末，一起剥皮案震惊了整个滨河市，随后出现的几起案子，更是在滨河造成了极大的恐慌，老刑警刘岩，带你破解...
    沈念sama阅读 212,884评论 6赞 492
  • 死咒
    序言：滨河连续发生了三起死亡事件，死亡现场离奇诡异，居然都是意外死亡，警方通过查阅死者的电脑和手机，发现死者居然都...
    沈念sama阅读 90,755评论 3赞 385
  • 救了他两次的神仙让他今天三更去死
    文/潘晓璐 我一进店门，熙熙楼的掌柜王于贵愁眉苦脸地迎上来，“玉大人，你说我怎么就摊上这事。” “怎么了？”我有些...
    开封第一讲书人阅读 158,369评论 0赞 348
  • 道士缉凶录：失踪的卖姜人
    文/不坏的土叔 我叫张陵，是天一观的道长。 经常有香客问我，道长，这世上最难降的妖魔是什么？ 我笑而不...
    开封第一讲书人阅读 56,799评论 1赞 285
  • 港岛之恋（遗憾婚礼）
    正文 为了忘掉前任，我火速办了婚礼，结果婚礼上，老公的妹妹穿的比我还像新娘。我一直安慰自己，他们只是感情好，可当我...
    茶点故事阅读 65,910评论 6赞 386
  • 恶毒庶女顶嫁案：这布局不是一般人想出来的
    文/花漫 我一把揭开白布。 她就那样静静地躺着，像睡着了一般。 火红的嫁衣衬着肌肤如雪。 梳的纹丝不乱的头发上，一...
    开封第一讲书人阅读 50,096评论 1赞 291
  • 城市分裂传说
    那天，我揣着相机与录音，去河边找鬼。 笑死，一个胖子当着我的面吹牛，可吹牛的内容都是我干的。 我是一名探鬼主播，决...
    沈念sama阅读 39,159评论 3赞 411
  • 双鸳鸯连环套：你想象不到人心有多黑
    文/苍兰香墨 我猛地睁开眼，长吁一口气：“原来是场噩梦啊……” “哼！你这毒妇竟也来了？” 一声冷哼从身侧响起，我...
    开封第一讲书人阅读 37,917评论 0赞 268
  • 万荣杀人案实录
    序言：老挝万荣一对情侣失踪，失踪者是张志新（化名）和其女友刘颖，没想到半个月后，有当地人在树林里发现了一具尸体，经...
    沈念sama阅读 44,360评论 1赞 303
  • 护林员之死
    正文 独居荒郊野岭守林人离奇死亡，尸身上长有42处带血的脓包…… 初始之章·张勋 以下内容为张勋视角 年9月15日...
    茶点故事阅读 36,673评论 2赞 327
  • 白月光启示录
    正文 我和宋清朗相恋三年，在试婚纱的时候发现自己被绿了。 大学时的朋友给我发了我未婚夫和他白月光在一起吃饭的照片。...
    茶点故事阅读 38,814评论 1赞 341
  • 活死人
    序言：一个原本活蹦乱跳的男人离奇死亡，死状恐怖，灵堂内的尸体忽然破棺而出，到底是诈尸还是另有隐情，我是刑警宁泽，带...
    沈念sama阅读 34,509评论 4赞 334
  • 日本核电站爆炸内幕
    正文 年R本政府宣布，位于F岛的核电站，受9级特大地震影响，放射性物质发生泄漏。R本人自食恶果不足惜，却给世界环境...
    茶点故事阅读 40,156评论 3赞 317
  • 男人毒药：我在死后第九天来索命
    文/蒙蒙 一、第九天 我趴在偏房一处隐蔽的房顶上张望。 院中可真热闹，春花似锦、人声如沸。这庄子的主人今日做“春日...
    开封第一讲书人阅读 30,882评论 0赞 21
  • 一桩弑父案，背后竟有这般阴谋
    文/苍兰香墨 我抬头看了看天上的太阳。三九已至，却和暖如春，着一层夹袄步出监牢的瞬间，已是汗流浃背。 一阵脚步声响...
    开封第一讲书人阅读 32,123评论 1赞 267
  • 情欲美人皮
    我被黑心中介骗来泰国打工， 没想到刚下飞机就差点儿被人妖公主榨干…… 1. 我叫王不留，地道东北人。 一个月前我还...
    沈念sama阅读 46,641评论 2赞 362
  • 代替公主和亲
    正文 我出身青楼，却偏偏与公主长得像，于是被迫代替她去往敌国和亲。 传闻我的和亲对象是个残疾皇子，可洞房花烛夜当晚...
    茶点故事阅读 43,728评论 2赞 351

推荐阅读更多精彩内容