数值分析第6章线性代数方程组的直接法.ppt

第六章,线性方程组 的直接解法,问题驱动：投入产出分析,投入产出分析是20世纪30年代由美国经济学家首先提出的，它是研究整个经济系统各部门之间“投入”与“产出”关系的线性模型，一般称为投入产出模型。国民经济各个部门之间存在着相互依存的关系，每个部门在运转中将其它部门的成品或半成品经过加工（称为投入）变为自己的产品（称为产出），如何根据各部门之间的投入-产出关系，确定各部门的产出水平，以满足社会的需求，是投入产出综合平衡模型研究的问题，试讨论如下简化问题。,设国民经济仅由农业、制造业和服务业三个部门构成，已知某年它们之间的投入和产出关系、外部需求、初始投入等如表所示（数字表示产值，单位为亿元）。,表6.1.1 国民经济各个部门间的关系,表中第一行数字表示农业总产出为100亿元，其中15亿元农产品用于农业生产本身，20亿元用于制造业，30亿元用于服务业，剩下的35亿元农产品用于满足外部需求。类似地可以解释第二、三行数字。第一列数字中，15亿元如前所述，30亿元是制造业对农业的投入，20亿元是服务业对农业的投入，35亿元的初始投入包括工资、税收、进口等，总投入100亿元和总产出相等。假定每个部门的产出和投入是成正比的，由表能够确定这三个部门的投入产出表，如表所示。,表6.1.2 投入产出表,表中的第一行，第二列的数字表示生产1个单位产值的制造业产品需要投入0.10个单位的产值的农产品，同样第三行、第一列的数字表示，生产1个单位产值的农产品需要0.20个单位的服务业产值。表的数字称为投入系数和消耗系数，如果技术水平没有变化，可以假设投入系数是常数。已知投入系数如表所示，若今年对农业、制造业和服务业的外部需求分别为50、150、100亿元，试计算三个部门的总产出分别为多少？,若共有n个部门，记一定时期内第i个部门的总产出为 xi,其中对第 j个部门的投入为xij，满足的外部需求为 di，则,(6.1.1),记第j个部门的单位产出需要第i个部门的投入为aij，在每个部门的产出与投入成正比的假定下，有,(6.1.2),投入系数即为aij，将(6.1.2)式代入(6.1.1)）式得方程组,用矩阵表示为,因此投入产出模型最终可归结为求解线性方程组的问题，下面介绍求解线性方程组数值方法。,AX=b,（3.1）,线性方程组数值解法的分类,直接法（适用于中等规模的n阶线性方程组）Gauss消去法及其变形 矩阵的三角分解法,迭代法（适用于高阶线性方程组）Jacobi迭代法 Gauss-Seidel迭代法 逐次超松弛法 共轭斜量法,1 高斯消去法,1三角形方程组的解法-回代法,（3.2）,（3.3）,2顺序高斯消去法,基本思想：通过消元将上述方程组 化为三角形方程组进行求解。,消元公式,回代公式,顺序Gauss消去法可执行的前提,定理 1 给定线性方程组，如果n阶方阵 的所有顺序主子式都不为零，即 则按顺序Gauss消去法所形成的各主元素 均不为零，从而Gauss 消去法可顺利执行。,注：当线性方程组的系数矩阵为对称正定或严格对角占优阵时，按Gauss消去法计算是稳定的。,3、列主元Gauss消去法计算步骤：,1、输入矩阵阶数n，增广矩阵 A(n,n+1);,2、对于,(1)按列选主元：选取 l 使,(2)如果，交换 A(n,n+1)的第k行与第l 行元素,(3)消元计算：,3、回代计算,4无回代过程的主元消去法(Gauss-Jordan),第一步：选主元，在第一列中选绝对值最大的元素，设第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消元法实际上就是将 方程组的系数矩阵化为行最简形矩阵。,5无回代消去法的应用,（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消元过程,。,2 解三对角方程组的追赶法,3 矩阵的三角分解法,高斯消元法的矩阵形式,每一步消去过程相当于左乘初等下三角矩阵Lk,记,A 的 LU 分解(LU factorization),定理2：（矩阵的三角分解）设A为n n实矩阵，如果 求解AX=b用顺序高斯消去法能够完成（即），则矩阵A可分解 为单位下三角矩阵L与上三角矩阵U的乘积。A=LU 且这种分解是唯一的。,注：（1）L 为单位下三角阵而 U 为一般上三角阵的分解称为Doolittle 分解（2）L 为一般下三角阵而 U 为单位上三角阵的分解称为Crout 分解。,Doolittle分解法：,LU 分解求解线性方程组,直接三角分解法解AX=b的计算公式,对于r=2,3,n计算,（2）计算U的第r行元素,（3）计算L的第r 列元素(r n),(1),(4),(5),例 用矩阵的三角分解法解方程组,Doolittle分解法的变形,紧凑格式的Doolittle分解法,例,所以,紧凑格式的列主元Doolittle分解法,例,4 平方根法,1矩阵的LDR分解,定理3：如果n阶矩阵A的所有顺序主子式均不等于零，则矩阵A存在唯一的分解式A=LDR，其中L和R分别是n阶单位下三角阵和单位上三角阵，D是对角元素不为零的n阶对角阵，上述分解称为A的LDR分解。,2平方根法,如果A为对称正定矩阵，则存在一个实的非奇异下三角矩阵，使A=LLT，且当限定的对角元素为正时，这种分解是唯一的，称为矩阵A的cholesky分解。,定理4：（对称正定矩阵的三角分解）,将对称 正定阵 A 做 LU 分解,即,则 仍是下三角阵，且有,对称正定阵cholesky分解的实现,用平方根法解对称正定线性代数方程组的算法,（1）对矩阵A进行Cholesky分解，即A=LLT，由矩阵乘法：,对于 i=1,2,n 计算,（2）求解下三角形方程组,（3）求解LTX=y,3改进平方根法,其中,改进平方根法解对称正定方程组的算法,令LTX=y，先解下三角形方程组 LDY=b 得,解上三角形方程组 LTX=Y 得,求解 时，A 和 的误差对解 有何影响？,设 A 精确，有误差，得到的解为，即,绝对误差放大因子,又,相对误差放大因子,5 线性方程组的性态和解的误差分析,6 Error Analysis for.,设 精确，A有误差，得到的解为，即,Wait a minute Who said that(I+A1 A)is invertible?,(只要 A充分小，使得,大,cond(H2)=,27,cond(H3),748,cond(H6)=,2.9 106,cond(Hn)as n,注：现在用Matlab数学软件可以很方便求矩阵的状态数!,定义2:设线性方程组的系数矩阵是非奇异的,如果cond(A)越大,就称这个方程组越病态.反之,cond(A)越小,就称这个方程组越良态.,一般判断矩阵是否病态，并不计算A1，而由经验得出。行列式很大或很小（如某些行、列近似相关）；元素间相差大数量级，且无规则；主元消去过程中出现小主元；特征值相差大数量级。,近似解的误差估计及改善：,设 的近似解为，则一般有,cond(A),改善方法(1)：,Step 1:近似解,Step 2:,Step 3:,Step 4:,若 可被精确解出，则有 就是精确解了。,经验表明：若 A 不是非常病态（例如：），则如此迭代可达到机器精度；但若 A 病态，则此算法也不能改进。,改善方法(2),对方程组进行预处理，即适当选择非奇异对角阵D，C,使求解 Ax=b 的问题转化为求解等价方程组 DACC-1x=Db,且使DAC 的条件数得到改善。（P88，例3.10）,用双精度进行计算，以便改善和减轻病态矩阵的影响。,历史与注记,高斯(Carl Friedrich Gauss,1777-1855)高斯是德国数学家、物理学家、天文学家。1777年生于德国布伦瑞克，1855年在哥廷根逝世。高斯是近代数学奠基者之一，在历史上影响之大，可以和阿基米德、牛顿、欧拉并列，有“数学王子”之称。1795 年进入格丁根大学学习。1798 年转入黑尔姆施泰特大学，1799 年获博士学位。,1807 年以后一直在格丁根大学任教授和哥廷根天文台台长，一直到逝世。1833年和物理学家W.E.韦伯共同建立地磁观测台，组织磁学学会以联系全世界的地磁台站网。,高斯的数学研究几乎遍及所有领域，在数论、代数学、非欧几何、复变函数和微分几何等方面都做出了开创性的贡献。高斯还把数学应用于天文学、大地测量学和磁学的研究，他在数值计算方面的贡献主要有：提出了高斯消去法和高斯迭代法、高斯积分，发明了最小二乘法原理等。,1810年，高斯提出了消去法，矩阵是凯莱在1855年提出的，而将高斯消去表示成矩阵分解是在20世纪40年代由Dwyer、冯.诺伊曼等人提出的。,参考文献,1 J.H.Wilkinson.Error analysis of direct methods of matrix inversion.J.ACM,8:281-330,1961.2 J.H.Wilkinson.Rounding Errors in Algebraic Processs.Prentice Hall,Englewood Cliffs,NJ,1963.3 J.H.Wilkinson.The Algebraic Eigenvalue Problem.Oxford University Press,New York,1965.4 G.E.Forsythe and C.B.Moler.Computer Solution of Linear Algebraic Systems.Prentice Hall,Englewood Cliffs,NJ,1967.,