第七章 线性代数方程组的直接法

第七章线性方程组的直接解法问题驱动：投入产出分析

投入产出分析是20世纪30年代由美国经济学家首先提出的，它是研究整个经济系统各部门之间“投入”与“产出”关系的线性模型，一般称为投入产出模型。国民经济各个部门之间存在着相互依存的关系，每个部门在运转中将其它部门的成品或半成品经过加工（称为投入）变为自己的产品（称为产出），如何根据各部门之间的投入-产出关系，确定各部门的产出水平，以满足社会的需求，是投入产出综合平衡模型研究的问题，试讨论如下简化问题。

设国民经济仅由农业、制造业和服务业三个部门构成，已知某年它们之间的投入和产出关系、外部需求、初始投入等如表7.1.1所示（数字表示产值，单位为亿元）。表7.1.1国民经济各个部门间的关系表中第一行数字表示农业总产出为100亿元，其中15亿元农产品用于农业生产本身，20亿元用于制造业，30亿元用于服务业，剩下的35亿元农产品用于满足外部需求。类似地可以解释第二、三行数字。第一列数字中，15亿元如前所述，30亿元是制造业对农业的投入，20亿元是服务业对农业的投入，35亿元的初始投入包括工资、税收、进口等，总投入100亿元和总产出相等。假定每个部门的产出和投入是成正比的，由表7.1.1能够确定这三个部门的投入产出表，如表7.1.2所示。表7.1.2投入产出表表中的第一行，第二列的数字表示生产1个单位产值的制造业产品需要投入0.10个单位的产值的农产品，同样第三行、第一列的数字表示，生产1个单位产值的农产品需要0.20个单位的服务业产值。表7.1.2的数字称为投入系数和消耗系数，如果技术水平没有变化，可以假设投入系数是常数。已知投入系数如表7.1.2所示，若今年对农业、制造业和服务业的外部需求分别为50、150、100亿元，试计算三个部门的总产出分别为多少？若共有n个部门，记一定时期内第i个部门的总产出为xi,其中对第j个部门的投入为xij，满足的外部需求为di，则

(1)记第j个部门的单位产出需要第i个部门的投入为aij，在每个部门的产出与投入成正比的假定下，有(2)投入系数即为aij，将(2)式代入(1)）式得方程组

用矩阵表示为或(3)因此投入产出模型最终可归结为求解线性方程组的问题，下面介绍求解线性方程组数值方法。AX=b（7.1）

线性方程组数值解法的分类直接法（适用于中等规模的n阶线性方程组）

◆Gauss消去法及其变形

◆矩阵的三角分解法迭代法（适用于高阶线性方程组）

◆Jacobi迭代法

◆

Gauss-Seidel迭代法

◆逐次超松弛法

◆共轭斜量法§1高斯消去法1．三角形方程组的解法---回代法（7.2）（7.3）

2．顺序高斯消去法

基本思想：通过消元将上述方程组化为三角形方程组进行求解。解n阶方程组的顺序高斯消去法的一般步骤：Step1：第一次消元消元因子Stepk：第k次消元消元因子共进行？步n

1三角方程组消元过程回代过程求解得：主元素注：若第一个主元素为0，因A为非奇异矩阵，可把第一列的非零元素所在的行调整到第一行即可。高斯消去法仍然可用。消元公式回代公式顺序Gauss消去法可执行的前提引理给定线性方程组，按顺序Gauss消去法所形成的各主元素均不为零，从而Gauss消去法可顺利执行的充要条件是n阶方阵的所有顺序主子式都不为零，即推论：如果A的顺序主子式不等于0，则(k=2,3,…,n)定理7.2：若A的所有顺序主子式

/determinantofleadingprincipalsubmatrices/均不为0，则高斯消元无需换行即可进行到底，得到唯一解。定理7.1：如果A非奇异，即A1存在，则可通过逐次消元及行交换，将方程组化为三角形方程组，求出唯一解。顺序Gauss消去法的计算量消元过程：k乘除法加减法1(n-1)+n(n-1)n(n-1)2(n-2)+(n-1)(n-2)(n-1)(n-2)﹕﹕﹕n-11+22总计其中，乘除法计算量总计：加减法计算量总计：回代过程：乘除法计算量：加减法计算量：引入主元素法的原因：举例用高斯消去法解方程组用四位浮点数进行计算，精确解舍入到4位有效数字为3.高斯主元素消去法解：［方法１］用顺序高斯消去法求解计算解为与精确解相比显然计算解是一个很坏的结果，不能作为方程组的近似解，其原因是：我们在消元计算时用了小主元0.001，使得约化后的方程组元素的数量级大大增长，使得在计算中发生严重的舍入误差，因此产生了较大的误差！［方法2］交换行，避免绝对值小的主元素作除数改进措施：对一般矩阵，最好每一步选取系数矩阵（或消元后的低阶矩阵）中绝对值最大的元素作为主元素，以使高斯消去法具有较好的数值稳定性！得计算解为本例启发：在采用高斯消去法解方程组时，小主元可能产生麻烦，故应避免采用绝对值小的主元素。完全主元素法设增广矩阵为第一步：首先在A中选取绝对值最大的元素作为主元素；然后交换到第一行、第一列的位置；再进行第一次消元，得矩阵(A(1)|b(1))→(A(2)|b(2))(1)消元过程第k步：在矩阵A(k)的右下方(n-k+1)阶子矩阵中选取绝对值最大的元素作为主元素；并通过行与列的互换将它换到第k行第k列的位置，然后进行第k次消元，得矩阵(A(k)|b(k))→(A(k+1)|b(k+1))第n-1步：经过n-1次消元，将原方程组化为其中,y1,y2,…,yn为未知数x1,x2,…,xn调换后的次序。(2)回代过程完全主元素消去法的特点：1.在选主元素时要花费较多机器时间。2.实时纪录x顺序的变化情况。3.实质为选主元的高斯消去法。

列主元素消去法选主元时仅考虑按列选取，然后换行使之变到主元位置上，再进行消元计算。列主元消去法的特点：（1）能够得到较高精度要求的解，数值稳定性好；（2）计算量较完全主元素法大大减少。列主元Gauss消去法计算步骤：1、输入矩阵阶数n，增广矩阵

A(n,n+1);2、对于(1)按列选主元：选取l使

(2)如果，交换A(n,n+1)的第k行与第l

行元素(3)

消元计算：3、回代计算高斯—约当消去法(Gauss-Jordan)基本思想：同时消去对角线上方和下方的元素。特点：(1)消元和回代同时进行；(2)乘除法的次数要比高斯消去法大。第一次消元(与高斯消去法不相同)第二次消元具体实现…故方程组的解为第一步：选主元，在第一列中选绝对值最大的元素，设第k行为主元行，将主元行换至第一行，将第一个方程中x1的系数变为1，并从其余n–1个方程中消去x1。第二步：在第二列后n–1个元素中选主元，将第二个方程中

x2的系数变为1，并从其它n–1个方程中消去x2。第k步：在第k列后n–k个元素中选主元，换行，将第k个方程xk的系数变为1，从其它n-1个方程中消去变量xk…………高斯—约当消去法步骤消元公式为：对k=1,2,…,按上述步骤进行到第n步后，方程组变为：即为所求的解注：Gauss-Jordan消元法实际上就是将方程组的系数矩阵化为行最简形矩阵。Gauss-Jordan消去法的应用（1）解线性方程组系设要解的线性方程组系为：AX=b1,AX=b2,…AX=bm上述方程组系可以写为AX=B=(b1,…,bm)因此 X=A-1B即为线性方程组系的解。

在计算机上只需要增加几组右端常数项的存贮单元，其结构和解一个方程组时一样。行系数右端（2）求逆矩阵设A=(aij)nn是非奇矩阵，A

0，且令由于AA-1=AX=I因此，求A-1的问题相当于解下列线性方程组相当于（1）中m=n,

B=I的情形。

（3）求行列式的值用高斯消去法将

A化成上三角形例

用Gauss-Jordan消去法解方程组

，并求出

其中解：把系数矩阵、单位矩阵和右端项组成增广矩阵，对增广矩阵实行Gauss-Jordan消元过程。故，系数矩阵的逆为矩阵的三角分解法高斯消去法的矩阵形式

每一步消去过程相当于左乘初等下三角矩阵Lk记其中，记其中，一般的，对A

的

LU

分解(LUfactorization)记则其中，于是矩阵三角分解的定义定义1设A为nn实矩阵，如果存在下三角矩阵L与上三角矩阵U，使得A=LU，则称之为矩阵A的三角分解；若存在单位下三角矩阵L，对角矩阵D及单位上三角矩阵R，使得A=LDR,则称之为矩阵A的LDR分解。注：（1）L为单位下三角阵而U为一般上三角阵的分解称为Doolittle

分解；（2）L为一般下三角阵而U为单位上三角阵的分解称为Crout分解。

定理2：设A为nn实矩阵，如果求解AX=b用顺序高斯消去法能够完成（即），则矩阵A可分解为单位下三角矩阵L与上三角矩阵U的乘积，即

A=LU且这种分解是唯一的。

矩阵三角分解的存在性直接三角分解法直接从矩阵A的元素得到计算L、U元素的递推公式，而不需任何中间步骤，这就是所谓的直接三角分解法。进行矩阵三角分解的意义若A实现了LU分解，则Ax=b(LU)x=bLy=bUx=y矩阵的三角分解法举例解方程组解：系数矩阵为于是求解原方程组可转化为如下两个三角形方程组：解得解得Ly=bUx=y???思考:能否不通过高斯消去法而直接获得系数矩阵的三角分解呢？即单位下三角阵上三角阵解：令由矩阵相等的定义得1矩阵的这种分解称为Doolittle分解杜利特尔

Doolittle分解法：通过比较法直接导出L和U的计算公式。思路LU分解求解线性方程组直接三角分解法解AX=b的计算公式对于r=2,3,…,n计算（2）计算U的第r行元素

（3）计算L的第r列元素(r

n)(1)(4)(5)例用矩阵的三角分解法解方程组Doolittle分解法的变形紧凑格式的Doolittle分解法例所以紧凑格式的列主元Doolittle三角分解法例Gauss消去法的变形1．矩阵的LDR分解定理3：如果n阶矩阵A的所有顺序主子式均不等于零，则矩阵A存在唯一的分解式A=LDR，其中L和R分别是n阶单位下三角阵和单位上三角阵，D是对角元素不为零的n阶对角阵，上述分解称为A的LDR分解。推论1(对称阵的三角分解定理Th7.7)设A为n阶对称阵，且A的所有顺序主子式均不为零，则A可唯一分解为A=LDLT其中，L为单位下三角阵，D为对角阵.证：因为A的所有顺序主子式不为零,故A有唯一的LU分解。将U再分解为其中，D为对角阵，U0为单位上三角阵，于是又由分解的唯一性得即对称正定阵定义及性质定义A为对称正定阵(3)A的特征值。(1)A是非奇异矩阵，且A-1亦是对称正定阵；(2)A的顺序主子式都大于零，即

1、性质2、对称正定阵的判别方法(充分条件)A为对称矩阵A为对称矩阵2．平方根法

如果A为对称正定矩阵，则存在一个实的非奇异下三角矩阵，使A=LLT，且当限定的对角元素为正时，这种分解是唯一的，称为矩阵A的cholesky分解。定理4：（对称正定矩阵的三角分解Th7.8）证明：由推论1可知A可分解为则由于对于任意非零向量x，y=(LT)-1x也为非零向量，于是由A的正定性即D正定，所以D的对角元素均为正数，令则将对称

正定阵

A做LU

分解U=uij=u11uij/uii111u22unn记为

A对称即记D1/2=则仍是下三角阵，且有对称正定阵cholesky分解的实现平方根法思路(与直接三角分解法解线性方程组类似推导公式)（1）将对称正定矩阵A分解为比较A与LLT的相应元素，可得计算公式分解唯一Th7.8设n阶对称正定矩阵A有分解

,先用待定系数法求L的元素用平方根法解线性代数方程组的计算步骤（1）对矩阵A进行Cholesky分解A=LLT，由矩阵乘法得：当j=1时，对j=2,3,…,n，由可得可得由（2）求解Ly=b

（3）求解LTx=yAx=b(LLT)x=bLy=bLTx=y用平方根法求解对称正定方程组时不需选取主元由可知因此因此，平方根法是数值稳定的。平方根法的稳定性和计算量平方根法计算量:n3/6（是LU分解的计算量的一半）。中间量得以控制，误差不会放大。优点：1、数值稳定。2、计算量小，大约为次乘除法，是一般矩阵A的LU分解计算量的一半。缺点：平方根法的优缺点：计算时要开平方。

设n阶对称正定矩阵A有分解

。其中L为单位下三角阵，D为对角阵。即

3．改进平方根法

=1

求解对称正定方程组Ax=b的改进平方根法(计算公式)：1．分解计算令为tij故有1、计算量小，同平方根方法一样，是一般矩阵A的

LU分解(消元法)计算量的一半，是目前解对称正定矩阵方程组的有效方法。3、精度较高。2、计算简单(没有开方运算)。改进平方根法的优点例用改进的平方根法求解方程组分析：显然系数矩阵为对称阵,再由TH7.7可进行LDLT分解解1）分解2）求解计算方程组的解：解三对角方程组的追赶法三对角方程组（1）对角占优矩阵:记三对角方程组（1）为:其中，追赶法设A为满足对角占优条件的n阶三对角阵，A=LU

—杜里特尔分解A=LU—Cront分解A=LDLT(一般三角分解)(追赶分解法)(平方根法)则有唯一三角分解A=LU由矩阵乘法，得:追赶法计算公式（1）分解计算公式()：（2）求解逆推公式

（3）求解逆推公式追赶法的本质追赶法的计算量和数值稳定性

追赶法的计算量比较小

追赶法是数值稳定的因为α,β,γ有界§5向量和矩阵的范数

一．向量范数定义1：设XRn，Ｘ

表示定义在Rn上的一个实值函数，称之为X的范数，它具有下列性质：（3）三角不等式：即对任意两个向量X、YRn，恒有

(1)非负性：即对一切XRn，X

0,Ｘ>0(2)齐次性：即对任何实数aR，XRn，在Rn上的向量x=(x1,…,xn)T∈Rn，三种常用的范数为：称为∞－范数或最大范数称为1－范数称为2－范数称为p－范数常用的向量范数举例：计算向量x=(1,-2,3)T的各种范数.解：向量范数的性质定义2：如果Rn中有两个范数||x||s与||x||t，存在常数m,M>0，使对任意n维向量x，有则称这两个范数等价.注：对两种等价范数而言，某向量序列在其中一种范数意义下收敛时，则在另一种范数意义下也收敛。定理1：定义在Rn上的向量范数

是变量X分量的

一致连续函数。定理2：在Rn上定义的任一向量范数

都与范数

等价，

即存在正数

M与m(M>m)

对一切XRn，不等式成立。推论：Rn上定义的任何两个范数都是等价的。

对常用范数，容易验证下列不等式：

定义3：设给定Rn中的向量序列{}，即其中若对任何i(i=1,2,…,n)都有则向量

称为向量序列{}的极限，或者说向量序列{}依坐标收敛于向量，记为定理3：向量序列{Xk}依坐标收敛于X*的充要条件是向量序列依范数收敛与依坐标收敛是等价的。二、矩阵范数1.矩阵范数的定义设对任意矩阵A∈Rn×n，按一定的规则有一实数与之对应，记为‖A‖，若‖A‖满足则称‖A‖为矩阵A的范数。正定性三角不等式矩阵范数与向量范数的相容性相容性定理4：设n

阶方阵A=(aij)nn，则（Ⅰ）与相容的矩阵范数是（Ⅱ）与相容的矩阵范数是其中1为矩阵ATA的最大特征值。（Ⅲ）与相容的矩阵范数是上述三种范数分别称为矩阵的1-范数、2-范数和∞-范数。2.常用的矩阵范数可以证明,对方阵和，有

(向量||·||2的直接推广)Frobenius范数:注：（1）（2）矩阵的Frobenius范数不是算子范数。举例：计算A的各种范数.解：下面计算2-范数令即故最大的特征值为所以例5:设A＝(aij)∈M.定义证明：这样定义的非负实数不是相容的矩阵范数.证明：设从而3．矩阵的范数与特征值之间的关系定理5：(Th7.16)矩阵A谱半径不超过A的任一矩阵范数，即

定义4：矩阵A的所有特征值的最大绝对值称为A的谱半径，记为：(Th7.17)如果A为对称矩阵，则

注:Rn×n中的任意两个矩阵范数也是等价的。定义6：设给定Rn×n中的矩阵序列{

}，若则称矩阵序列{}收敛于矩阵A，记为称为A与B之间的距离。定义5：设为上的矩阵范数，定理6（TH7.18）如果，则为非奇异矩阵，且其中，是指矩阵的算子范数。定理7（TH8.2）设B∈Rn×n，则由B的各幂次得到的矩阵序列Bk,k=0,1,2…)收敛于零矩阵，即的充要条件为

4.矩阵的条件数定义7设矩阵为非奇异矩阵，则称为矩阵

的条件数，其中是矩阵的算子范数。对矩阵的任意一个算子范数有(2)cond(kA)=cond(A),k为非零常数;(3)若，则注:

cond(A)与所取的范数有关常用条件数有：cond(A)2特别地，若A对称，则cond(A)1=‖A‖1‖‖1cond(A)=‖A‖

‖‖求解时，A和的误差对解有何影响？设A精确，有误差，得到的解为，即绝对误差放大因子又相对误差放大因子§6线性方程组的性态和解的误差分析§6ErrorAnalysisfor.

设精确，A有误差，得到的解为，即

Waitaminute…

Whosaidthat(I+A1A)isinvertible?(只要A充分小，使得是关键的误差放大因子，称为A的条件数，记为cond(A),越则A越病态，难得准确解。大例：Hilbert阵cond(H2)=27cond(H3)748cond(H6)=2.9106cond(Hn)asn注：现在用Matlab数学软件可以很方便求矩阵的状态数!定义8:设线性方程组的系数矩阵是非奇异的,如果cond(A)越大,就称这个方程组越病态.反之,cond(A)越小,就称这个方程组越良态.一般判断矩阵是否病态，并不计算A1，而由经验得出。行列式很大或很小（如某些行、列近似相关）；元素间相差大数量级，且无规则；主元消去过程中出现小主元；特征值相差大数量级。近似解的误差估计及改善：设的近似解为，则一般有cond(A)误差上限改善方法(1)：Step1:近似解Step2:Step3:Step4:若可被精确解出，则有就是精确解了。经验表明：若A不是非常病态（例如：），则如此迭代可达到机器精度；但若A病态，则此算法也不能改进。改善方法(2)对方程组进行预处理，即适当选择非奇异对角阵D，C,使求解Ax=b的问题转化为求解等价方程组