线性代数方程组的直接法.ppt
2023/2/201王王 淑淑 栋栋办办公室：公室：J13-409电话电话：88032786（H）E-mail：数 值 分 析（Numerical Analysis）2023/2/202第七章第七章 解线性方程组的解线性方程组的直接直接方法方法数 值 分 析（Numerical Analysis）AXAX=b b很多实际问题都可归结为求解线性方程组：很多实际问题都可归结为求解线性方程组：2023/2/203写成矩阵形式：写成矩阵形式：线性方程组数值解法的分类：线性方程组数值解法的分类：直接法：直接法：经过有限步算术运算即可求得方程组经过有限步算术运算即可求得方程组精精确解确解的方法的方法（适用于中等规模的（适用于中等规模的n阶线性方程组）阶线性方程组）Gauss消去法及其变形消去法及其变形矩阵的三角分解法矩阵的三角分解法 迭代法：迭代法：用某种极限过程去逐步逼近线性方程用某种极限过程去逐步逼近线性方程组组精确解精确解的方法。（的方法。（适用于高阶线性方程组适用于高阶线性方程组）Jacobi迭代法迭代法Gauss-Seidel迭代法迭代法逐次超松弛法逐次超松弛法2023/2/2042023/2/2051Gauss消去法消去法或或或或 AXAX=b b消元手续消元手续2023/2/206举例说明消去法的基本思想：举例说明消去法的基本思想：基本思想：用逐次消去未知数的方法把原来方程组基本思想：用逐次消去未知数的方法把原来方程组AXAX=b b化为化为与其等价的三角方程组，然后再用回代法写出方程组的解。事与其等价的三角方程组，然后再用回代法写出方程组的解。事实上，将方程组化为与其等价的三角方程组的过程就是用实上，将方程组化为与其等价的三角方程组的过程就是用行的行的初等变换初等变换将原来方程组的系数矩阵化为简单形式，然后再求解。将原来方程组的系数矩阵化为简单形式，然后再求解。1消元消元 基本思想：基本思想：通过消元手续将上述方程组化为三角形通过消元手续将上述方程组化为三角形方程组进行求解。方程组进行求解。2023/2/207Guass消去法消去法消元公式消元公式消元公式消元公式2023/2/2082.回代回代2023/2/209回代公式回代公式回代公式回代公式2023/2/2010Gauss消去法可执行的前提：消去法可执行的前提：【定理定理1】按顺序按顺序Gauss消去法所形成的各主元素消去法所形成的各主元素的充要条件是的充要条件是矩阵矩阵的的顺序主子式顺序主子式，即，即2023/2/2011【推论推论】如果如果的顺序主子式的顺序主子式2023/2/2012【定理定理2】如果如果的所有顺序主子式都不为零，即的所有顺序主子式都不为零，即则可通过则可通过Gauss消去法消去法将前面的方程组约化为将前面的方程组约化为三角三角方程组。方程组。2023/2/2013作业：作业：P1971，2，7不进行交换两行不进行交换两行 的初等变换！的初等变换！2023/2/2014矩阵的三角分解矩阵的三角分解借助矩阵理论分析消去法！建立借助矩阵理论分析消去法！建立Gauss消去法消去法和矩阵理论和矩阵理论间的关系！对于间的关系！对于或或或或 AXAX=b b在在Gauss消元法中，消元法中，每一步消去过程相当于每一步消去过程相当于左乘初等变换左乘初等变换矩阵矩阵Lk2023/2/20152023/2/2016A的的LU分解分解其中其中 单位下三角形。单位下三角形。2023/2/2017【定理定理3】（矩阵的（矩阵的LU分解）设分解）设A为为n n矩阵，如果矩阵，如果解解AX=b用用高斯消去法高斯消去法（限制不进行行的交换，即（限制不进行行的交换，即 ）能够完成，则矩阵）能够完成，则矩阵A可分解为可分解为单位下三角矩阵单位下三角矩阵L与与上三上三角矩阵角矩阵U的乘积，即的乘积，即 A=LU且这种分解是且这种分解是唯一唯一的。的。2023/2/2018注注注注：（1 1）L为单位下三角阵而为单位下三角阵而U为一般上三角为一般上三角阵的分解称为阵的分解称为Doolittle分解分解（2）L为一般下三角阵而为一般下三角阵而U为单位上三角为单位上三角阵的分解称为阵的分解称为Crout分解。分解。2023/2/20192023/2/2020系数矩阵系数矩阵 的三角分解！的三角分解！21计算量计算量（1）消元消元过程的计算量：第过程的计算量：第1步步计算乘数计算乘数mi1(i=1,2,n)，需要，需要n-1次除法运算；计算次除法运算；计算需要需要次乘法和次乘法和次加、减法。次加、减法。第第k步步 加、减法次数加、减法次数乘法次数乘法次数除法次数除法次数12n-1111合计合计22第第k步步 加、减法次数加、减法次数乘法次数乘法次数除法次数除法次数1232.21nn-1n-11合计合计(2)计算计算b(n)计算量计算量:进行进行次乘法和次乘法和加、减法；加、减法；(3)回代过程的计算量回代过程的计算量:总计算量：总计算量：次乘除法和次乘除法和次加减法次加减法较大时较大时乘除法次数：乘除法次数：较大时较大时加减法次数：加减法次数：2Gauss主元素主元素消去法消去法【例例4】用用Gauss消去法解方程组消去法解方程组选主元素的必要性！选主元素的必要性！要求用具有舍入的要求用具有舍入的10位浮点数进行计算。位浮点数进行计算。精确到精确到10位的精确解为：位的精确解为：【解法解法1】（高斯消去法高斯消去法）消元：消元：舍去或者舍去或者说被说被“吃吃”舍去或者说被舍去或者说被“吃吃”计算解：计算解：2023/2/2023显然显然,计算解与精确解解相差太大，原因是用很小的数计算解与精确解解相差太大，原因是用很小的数作除数，使得舍入误差太大，从而计算结果不可靠。作除数，使得舍入误差太大，从而计算结果不可靠。【解法解法2】用用行变换行变换的高斯消去法的高斯消去法消元：消元：计算解计算解:该该结结果果较较好好。该该例例子子说说明明，在在采采用用高高斯斯消消去去法法解解方方程程组组时时，应应避避免免采采用用绝绝对对值值很很小小主主元元素素。对对一一般般系系数数矩矩阵阵，最最好好保保持持乘乘数数，因因此此在在高高斯消去法中引进斯消去法中引进选主元素选主元素技巧。技巧。241、完全主元素消去法、完全主元素消去法选主元素消元法：选主元素消元法：为为非奇异矩阵非奇异矩阵，第一步：第一步：（1）选选主主元元素素：在在A中中选选取取绝绝对对值值最最大大的的元元素素作作为为主主元元素素，即即确定确定使使（2）交交换换行行列列：交交换换中中第第1行行和和i1行行元元素素，第第1列列和和第第j1列列元元素素。注注意意调调换换两两未未知知量量，并并做做记记录录。交交换换后后的的元元素素仍仍记为记为。（3）第第1次消元计算次消元计算：重复进行上述过程，设已完成第重复进行上述过程，设已完成第1步至第步至第k-1步步的选主元素、的选主元素、交换行列和消元计算，使交换行列和消元计算，使为增广矩阵。为增广矩阵。25第第k步：步：（1）选主元素：选取）选主元素：选取使使（2）交交换换行行列列：交交换换中中第第k行行和和ik行行元元素素，第第k列列和和第第jk列元素。列元素。（3）消元计算：）消元计算：2023/2/2026回代求解回代求解工作量大。工作量大。经过上述过程，方程组约化为：经过上述过程，方程组约化为：缺点：缺点：优点：数值稳定优点：数值稳定改进方法：改进方法：列列主元消去法，且此时主元消去法，且此时其中其中是未知数是未知数调换后的顺序，则调换后的顺序，则完全主元素消去法优缺点：完全主元素消去法优缺点：完全完全主元素消去法步骤见主元素消去法步骤见P172!27列列主主元元素素法法考考虑虑依依次次按按列列选选主主元元素素，然然后后换换行行，再再进进行行消消元元计计算算。设设已已完完成成第第1步步第第k-1步步计计算算，得得到到与与原原方方程组等价的方程组程组等价的方程组其中其中方框内为第方框内为第k步选步选主元主元素区域。素区域。2、列列主元素消去法主元素消去法28输入矩阵阶数输入矩阵阶数n，增广矩阵增广矩阵 A(n,n+1);对于对于(1)按列选主元：选取按列选主元：选取 l，使使(2)如果如果 ，交换，交换 A(n,n+1)的第的第k行与第行与第l行元素行元素(3)消元计算消元计算 ：回代计算回代计算2023/2/2029列主元素消去法步骤见列主元素消去法步骤见P173!矩阵运算矩阵运算描述列主元素消去法！描述列主元素消去法！30先换行再消元！先换行再消元！【定理定理4】（列主元素的三角分解）如果为非奇异矩阵，则存（列主元素的三角分解）如果为非奇异矩阵，则存在排列矩阵在排列矩阵P，使得，使得PA=LU，其中，其中L是是单位下单位下三角阵，三角阵，U是是上上三角三角阵。阵。【例例5】用用列列主元素消去法解方程组主元素消去法解方程组我们用我们用4位浮点数进行计算。位浮点数进行计算。解解：消元：消元：舍去或者说被舍去或者说被“吃吃”已知精确解：已知精确解：2023/2/20回代计算解：回代计算解：高斯高斯选主元消去法选主元消去法的步骤：的步骤：注：该解若取注：该解若取两位有效数字两位有效数字，则与精确解完全相同。，则与精确解完全相同。u优点：优点：数值稳定数值稳定u修正方法：修正方法：消元，回代。消元，回代。高斯高斯-约当（约当（Gauss-Jordam）消去法。消去法。u缺点：缺点：既消元，又回代既消元，又回代2023/2/20323、高斯、高斯约当（约当（GaussJordan)消去法消去法高高斯斯消消去去法法始始终终是是消消去去对对角角线线下下方方的的元元素素，GaussJordan消消去去法法要要消消去去对对角角线线下下方方和和上上方方的的元元素素。假假设设G-J消消去去法法已已完完成成第第1步步第第k-1步步，得得到到与与原原方方程组等价的方程组程组等价的方程组，其中，其中第第k步计算：步计算：（1）按列选主元：即确定）按列选主元：即确定ik，使使（2）换行：交换）换行：交换的第的第k行和第行和第ik行元素行元素（3）消元计算：）消元计算：33（4）计算主行（主元素所在行）计算主行（主元素所在行）计算解：计算解：上述过程结束后，有上述过程结束后，有2023/2/2034说说明明：在在解解方方程程组组，一一般般不不用用高高斯斯-约约当当消消去去法法。因因为为计计算算量量太太大大，但但是是在在解解多多个个方方程程组组而而它它们们的的系系数数矩矩阵阵相相同同时时，用用此此方方法，即求系数矩阵的逆矩阵法，即求系数矩阵的逆矩阵A-1时，比较合适。时，比较合适。设设为为非非奇奇异异矩矩阵阵。如如果果用用列列主主元元G-J消消去去法法将将（A,I）化化为（为（I,T），则，则。不用回代，将不用回代，将A化为化为单位矩阵单位矩阵，则解为常数项列。，则解为常数项列。【定理定理5】（高斯（高斯约当法求逆矩阵）约当法求逆矩阵）优点：优点：缺点：计算量较大，大约是缺点：计算量较大，大约是次乘除法。次乘除法。注注：该方法与高等代数中求逆矩阵方法的不同之处是有：该方法与高等代数中求逆矩阵方法的不同之处是有选主元选主元，实际上选主元就是交换两行的位置，仍是初等变换，在一般的实际上选主元就是交换两行的位置，仍是初等变换，在一般的求逆矩阵方法中也有交换两行元素。求逆矩阵方法中也有交换两行元素。2023/2/2035【例例6】用列主元素用列主元素G-J消去法求消去法求的逆矩阵的逆矩阵A-1解：解：2023/2/2036所以所以G-J列主元列主元求逆求逆算法见算法见P176算法算法3！2023/2/20372023/2/2038作业：作业：P199122023/2/20393Gauss消去法消去法的变形的变形1.直接直接三角分解法三角分解法直接从矩阵直接从矩阵A的元素得到计算的元素得到计算，U元素的递推公式，而不需元素的递推公式，而不需要任何中间步骤，称之为要任何中间步骤，称之为直接三角分解法直接三角分解法。这样，。这样，AX=bL(UX)=bLY=b，UX=Y。不选主元不选主元的三角分解的三角分解(Doolittle分解法分解法)通过比较法直接导出通过比较法直接导出L 和和 U 的计算公式。的计算公式。思思路路40通过通过n步可以直接计算定出步可以直接计算定出L，U的元素，其中第的元素，其中第r步确定步确定U的的第第r行和行和L的第的第r列元素。列元素。第第1步：步：设第设第r-1步已确定步已确定U的第的第r-1行和行和L的第的第r-1列元素。列元素。41第第r步：步：故：故：故：故：LU 分解求解线性方程组分解求解线性方程组:2023/2/2042直接三角分解法解直接三角分解法解AX=b的计算公式的计算公式:对于对于r=2,3,n，(2)计算计算U的第的第r行元素行元素(3)计算计算L的第的第r 列元素列元素(r n)(1)2023/2/2043(4)(5)2023/2/2044【例例7】用直接三角分解法解用直接三角分解法解解：用分解计算公式得：解：用分解计算公式得：求解：求解：2023/2/20452023/2/2046选主元的三角分解选主元的三角分解选主元素三角分解算法见选主元素三角分解算法见P179!2.平方根法平方根法(对称正定矩阵而言对称正定矩阵而言)矩阵的矩阵的LDR分解分解【定理定理6】若若n阶矩阵阶矩阵A的所有顺序主子式均不等于零，的所有顺序主子式均不等于零，则矩阵则矩阵A存在唯一的分解式存在唯一的分解式A=LDR，其中其中L和和R分别是分别是n阶阶单位单位下下三角阵和三角阵和单位单位上上三角阵，三角阵，D是是n阶对角元素阶对角元素不为零的对角阵，上述分解也称为不为零的对角阵，上述分解也称为A的的LDR分解分解。2023/2/2047平方根法平方根法如果如果A为为对称对称矩阵，矩阵，且且A的所有顺序主子式均不的所有顺序主子式均不等于零，等于零，则则A可唯一分解为可唯一分解为A=LDLT，其中其中L是是n阶阶单位单位下下三角阵和三角阵和D是是对角阵对角阵。【定理定理7】（对称对称矩阵的三角分解定理）矩阵的三角分解定理）2023/2/2048设设A是是对称对称 正定正定矩阵，矩阵，做做LU分解分解U=uij=u11uij/uii111u22unn记为记为 A 对称对称即即记记 D1/2=则则 仍是仍是下下三角阵，且有三角阵，且有2023/2/2049【定理定理5】（对称正定对称正定矩阵的三角分解或矩阵的三角分解或Cholesky分解）分解）2023/2/2050如果如果A为对称正定矩阵，则存在一个实的非奇异为对称正定矩阵，则存在一个实的非奇异下三下三角角矩阵，使矩阵，使A=LLT，且当限定的对角元素为正时，且当限定的对角元素为正时，这种分解是唯一的。这种分解是唯一的。如何确定如何确定L呢呢？下面用下面用直接法直接法确定确定L中元素。中元素。2023/2/2051其中其中，注意到，注意到所以所以因此因此用平方根法解线性代数方程组的算法用平方根法解线性代数方程组的算法（1）对矩阵对矩阵A进行进行Cholesky分解分解，即即A=LLT，由矩阵乘法由矩阵乘法：对于对于 i=1,2,n 计算计算2023/2/2052（2）求解下三角形方程组求解下三角形方程组 （3）求解求解LTX=y2023/2/2053改进平方根法改进平方根法 其中其中2023/2/2054改进平方根法改进平方根法解对称正定方程组的算法解对称正定方程组的算法 2023/2/2055令令LTX=y，先解下三角形方程组先解下三角形方程组LDY=b 得得解上三角形方程组解上三角形方程组 LTX=Y 得得 2023/2/20563.追赶法追赶法2023/2/20572023/2/20584向量和矩阵的范数向量和矩阵的范数 1向量向量范数范数【定义定义】设设X R n，X 表示定义在表示定义在Rn上的上的一个实值函数一个实值函数，称之为称之为X的范数的范数，它具有下列性质它具有下列性质：（3）三角不等式三角不等式：即对任意两个向量即对任意两个向量X,Y R n，恒有恒有(1)非负性非负性：即对一切即对一切X R n，X 0,X 0(2)齐次性齐次性：即对任何实数即对任何实数a R，X R n，2023/2/2059由上可得：由上可得：设设X=(x1,x2,xn)T，则有则有（1）（2）（3）三个常用的向量范数：三个常用的向量范数：范数等价范数等价:设设s 和和t是是Rn上任意两种范数，若存在上任意两种范数，若存在常数常数c1，c2 0使得使得:,则称则称s 和和t 等价。等价。2023/2/2060【定理定理5】(N(x)的连续性的连续性)定义在定义在Rn上的向量范数上的向量范数是是X分量分量的的连续函数。连续函数。【定理定理6】（向量范数的等价性）在（向量范数的等价性）在Rn上定义的任一向量范数上定义的任一向量范数都与范数都与范数等价，即存在正数等价，即存在正数M与与m(Mm)对一切对一切X Rn，不等式不等式成立。成立。推论：推论：Rn上定义的任何两个范数都是等价的。上定义的任何两个范数都是等价的。2023/2/2061对常用范数，容易验证下列不等式：对常用范数，容易验证下列不等式：2023/2/2062【定义定义2】设设 为为Rn中一向量序列中一向量序列，即即其中其中若对任何若对任何i(i=1,2,n)都有都有则向量则向量 称为向量序列称为向量序列 的极限的极限，或者说向量序列或者说向量序列依坐标收敛于向量依坐标收敛于向量 ，记为记为2023/2/2063【定理定理7.13】向量序列向量序列依坐标收敛于依坐标收敛于X*的充要条件是的充要条件是向量序列依范数收敛与依坐标收敛是等价的。向量序列依范数收敛与依坐标收敛是等价的。【定义定义3】设设A为为n 阶方阵，阶方阵，Rn中已定义了向量范数中已定义了向量范数，则称则称为矩阵为矩阵A的范数或模，的范数或模，记为记为。即。即2023/2/2064（1）当）当A=0时，时，0，当，当A 0时，时，0（2）对任意实数对任意实数k和任意和任意A，有，有（3）对任意两个对任意两个n阶矩阵阶矩阵A，B有有（4）对任意两个对任意两个n阶矩阵阶矩阵A，B，有有2023/2/20652矩阵的范数矩阵的范数【定义定义3】（矩阵范数）如果矩阵的（矩阵范数）如果矩阵的某个非负的实值某个非负的实值函数函数满足条件：满足条件：则称则称 是是 上的一个矩阵范数上的一个矩阵范数(模模)。2023/2/2066矩阵和向量相容的范数：矩阵和向量相容的范数：其中其中满足矩阵范数的定义，称为满足矩阵范数的定义，称为A的的Frobenius范数。范数。【例例5】设设A(aij)M。定义定义证明：这样定义的非负实数不是证明：这样定义的非负实数不是相容相容的矩阵范数的矩阵范数.证明：设证明：设从而从而2023/2/2067【定理定理8】设设n 阶方阵阶方阵A=(aij)n n，则，则（）与）与相容的矩阵范数是相容的矩阵范数是（）与）与相容的矩阵范数是相容的矩阵范数是其中其中 1为矩阵为矩阵ATA的最大特征值。的最大特征值。（）与）与相容的矩阵范数是相容的矩阵范数是上述三种范数分别称为矩阵的上述三种范数分别称为矩阵的1-范数，范数，2-范数和范数和-范数。范数。2023/2/2068可以证明可以证明,对方阵对方阵和和，有，有(向量向量|2的直接推广的直接推广)Frobenius范数范数:2023/2/20693矩阵矩阵范数与特征值之间的关系范数与特征值之间的关系【定理定理9】矩阵矩阵A 的任一特征值的绝对值不超过的任一特征值的绝对值不超过A的范数，即的范数，即 【定义定义4】矩阵矩阵A的诸特征值的诸特征值 的最大绝对值称为的最大绝对值称为A的的谱半径谱半径，记为：记为：并且如果并且如果A A为对称矩阵，则为对称矩阵，则 2023/2/2070注注:Rnn中的任意两个矩阵范数也是等价的。中的任意两个矩阵范数也是等价的。【定义定义5】设设|为为Rnn上的矩阵范数，上的矩阵范数，A,BRnn称称|A-B|为为A与与B之间的距离。之间的距离。【定义定义6】设给定设给定Rnn中的矩阵序列中的矩阵序列，若，若则称矩阵序列则称矩阵序列收敛于矩阵收敛于矩阵A，记为，记为2023/2/2071【定理定理10】设设BRnn，则由，则由B的各幂次得到的的各幂次得到的矩阵序列矩阵序列Bk,k=0,1,2)收敛于零矩阵收敛于零矩阵的充要条件为的充要条件为。2023/2/2072求解求解时，时，A和和的误差对解的误差对解有何影响？有何影响？设设A精确，精确，有误差有误差，得到的解为，得到的解为，即，即绝对误差放大因子绝对误差放大因子又又相对误差放大因子相对误差放大因子5 误差分析误差分析2023/2/2073 设设精确，精确，A有误差有误差，得到的解为，得到的解为，即，即 Wait a minute Who said that(I+A 1 A)is invertible?(只要只要 A充分小，使得充分小，使得是关键是关键的误差放大因子，称为的误差放大因子，称为A的条件数，记为的条件数，记为cond(A),越越则则A 越病态，越病态，难得准确解。难得准确解。大大2023/2/2074【定义定义7】设设A为为n 阶非奇异矩阵，称数阶非奇异矩阵，称数 为矩阵为矩阵A的的条件数。条件数。条件数的性质：条件数的性质：）cond(A)1）cond(kA)=cond(A),k 为非零常数为非零常数）若）若 ，则，则2023/2/2075注注:cond(A)与所取的范数有关与所取的范数有关常用条件数有：常用条件数有：cond(A)2特别地，若特别地，若A 对称，则对称，则cond(A)1=A1 1cond(A)=A 2023/2/2076【例例】Hilbert阵阵cond(H2)=27cond(H3)748cond(H6)=2.9 106cond(Hn)当当n 注：现在用注：现在用Matlab数学软件可以很方便求矩阵的状态数数学软件可以很方便求矩阵的状态数!【定义定义2】设线性方程组的系数矩阵是非奇异的，如果设线性方程组的系数矩阵是非奇异的，如果cond(A)越大，就称这个方程组越病态。反之，越大，就称这个方程组越病态。反之，cond(A)越小，越小，就称这个方程组越良态。就称这个方程组越良态。2023/2/2077一般判断矩阵是否病态，并不计算一般判断矩阵是否病态，并不计算A 1，而由经验得出。，而由经验得出。行列式很大或很小（如某些行、列近似相关）；行列式很大或很小（如某些行、列近似相关）；元素间相差大数量级，且无规则；元素间相差大数量级，且无规则；主元消去过程中出现小主元；主元消去过程中出现小主元；特征值相差大数量级。特征值相差大数量级。2023/2/2078 近似解的误差估计及改善：近似解的误差估计及改善：设设的近似解为的近似解为，则一般有，则一般有cond(A)误差上限误差上限改善方法改善方法(1)：Step 1:近似解近似解Step 2:Step 3:Step 4:若若可被精确解出，则有可被精确解出，则有就是精确解了。就是精确解了。经验表明经验表明：若：若A 不是非常病态（例如：不是非常病态（例如：），则如），则如此迭代可达到机器精度；但若此迭代可达到机器精度；但若A 病态，则此算法也不能改病态，则此算法也不能改进。进。2023/2/2079 改善方法改善方法(2)(2)对方程组进行预处理，即适当选择非奇异对角阵对方程组进行预处理，即适当选择非奇异对角阵D D，C,C,使求解使求解 Ax=b Ax=b 的问题转化为求解等价方程组的问题转化为求解等价方程组 DACCDACC-1-1x=Db,x=Db,且使且使DACDAC 的条件数得到改善。（的条件数得到改善。（P88P88，例，例3.103.10）用双精度进行计算，以便改善和减轻病态矩阵的影响用双精度进行计算，以便改善和减轻病态矩阵的影响。2023/2/2080
