数值分析03-线性方程组直接解法.ppt

阜师院数科院第三章 线性方程组直接解法,3-1,第三章,线性方程组直接解法,阜师院数科院第三章 线性方程组直接解法,3-2,第三章目录,1.Gauus 消元法2.主元素法 2.1 引入主元素法的必要性 2.2 列主元素法 2.3 全主元素法 2.4 解三对角方程组的追赶法3.矩阵分解法 3.1 Gauss消去法的矩阵形式 3.2 矩阵的三角分解 3.3 直接三角分解法4.平方根法与改进的平方根法5.矩阵求逆6.方程组的性态和条件数,阜师院数科院第三章 线性方程组直接解法,3-3,线性方程组的概念,在科学研究和工程技术中所提出的计算问题中，线性方程组的求解问题是基本的，常见的，很多问题如插值函数，最小二乘数据拟合，构造求解微分方程的差分格式等，都包含了解线性方程组问题，因此，线性方程组的解法在数值计算中占有较重要的地位。,设n阶线性方程组：,其矩阵形式为：,Ax=b(3-2),其中：,阜师院数科院第三章 线性方程组直接解法,3-4,求解Ax=b，曾经学过高斯（Gauss）消元法，克莱姆（Cramer）法则，矩阵变换法等，但已远远满足不了实际运算的需要，主要体现两个方面：一是运算的快速和准确，其次是方程组的个数增大时的计算问题。如何建立能在计算机上可以实现的有效而实用的解法，具有极其重要的意义，我们也曾指出过，Cramer法则在理论上是绝对正确的，但当n较大时，在实际计算中却不能用。,线性方程组的概念（续）,如果线性方程组Ax=b的系数行列式不为零，即det(A)0，则该方程组有唯一解。,阜师院数科院第三章 线性方程组直接解法,3-5,线性方程组的数值解法,解线性方程组的数值方法大致分为两类：,请注意：由于在计算中某些数据实际上只能用有限位小 数，即不可避免地存在着舍入误差的影响，因 而即使是确解法，也只能求到近似解。直接法在求解中小型线性方程组（100个），特别是系数矩阵为稠密型时，是常用的、非常好的方法。,直接法：指假设计算过程中不产生含入误差，经过有 限步四则运算可求得方程组准确解的方法。,2.迭代法：从给定的方程组的一个近似值出发，构造某种算法逐步将其准确化，一般不能在有限步内得到准确解。,这一章介绍计算机上常用的直接法，它们都是以Gauss消元法为基本方法，即先将线性方程组化为等价的三角形方程组，然后求解。,阜师院数科院第三章 线性方程组直接解法,3-6,1 Gauss消元法,Gauss消元法是最基本的一种方法，下例说明其基本思想:,例1,解线性方程组：,解：消去x1，进行第一次消元：首先找乘数，以-12乘第一个方程加到第二个方程，以18乘第一个方程加到第三个方程上可得同解方程组：,阜师院数科院第三章 线性方程组直接解法,3-7,例1（续）,上述Gauss消元法的基本思想是：先逐次消去变量，将方程组化成同解的上三角形方程组，此过程称为消元过程。然后按方程相反顺序求解上三角形方程组，得到原方程组的解，此过程称为回代过程。,再消一次元得：,二次消元后将方程化为倒三角形式，然后进行回代容易解出：x3=3,x2=2,x1=1。,我们的目的，是要总结归纳出一般情况下的n阶线性方程组的消元公式和回代求解公式，从而得到求解n阶线性方程组的能顺利在计算机上实现的行之有效的算法。,阜师院数科院第三章 线性方程组直接解法,3-8,Gauss消元法的基本步骤1（4阶）,为能更清楚地得到算法，下面以4阶线性方程组为例总结求解步骤，并且很容易地可推广至一般的n阶线性方程组。,阜师院数科院第三章 线性方程组直接解法,3-9,可以检查，分别以li1乘第一个方程加到第i个方程上可以完成第一次消元，得同解方程组：,变化以后的方程组系数及右边的常数项可总结出如下的计算公式：,Gauss消元法的基本步骤2（4阶）,阜师院数科院第三章 线性方程组直接解法,3-10,Gauss消元法的基本步骤3（4阶）,以方程组中第i个方程减去第二个方程乘li2(i=3,4)，完成第二次消元。,上标为3的系数和右端项可由下面公式计算：,阜师院数科院第三章 线性方程组直接解法,3-11,Gauss消元法的基本步骤4（4阶）,第三步：消元（4阶方程组需进行3次消元）将上述 A(3)X=b(3)中最后一个方程中的x3消为零：,然后可回代求解：由于A(4)为上三角形，所以可按变量的逆序逐步回代求原方程组的解：,上述 消元、回代求解过程很容易推广到一般的n阶线性方程组。,经过上述消元步骤，得到同解的上三角形方程组：A(4)x=b(4),阜师院数科院第三章 线性方程组直接解法,3-12,Gauss消元法的消元过程1、2（n阶）,一般地，设 n阶方程组：,消元过程为：,阜师院数科院第三章 线性方程组直接解法,3-13,Gauss消元法的消元过程3（n阶）,第k步消元后同解方程组中上标为k+1的元素的计算公式见下屏,阜师院数科院第三章 线性方程组直接解法,3-14,照此消元下去，完成n1次消元后，可将原方程组化成同解的上三角形方程组如下：,Gauss消元法的消元过程3（n阶）,阜师院数科院第三章 线性方程组直接解法,3-15,Gauss消元法的回代过程（n阶）,回代过程：逐步回代求得原方程组的解,阜师院数科院第三章 线性方程组直接解法,3-16,Gauss消元法的计算量,由于在计算机中作乘除运算量所需时间远大于作加减运算所需时间，故只考虑作乘除运算量。由消元法步骤知，第k次消元需作nk次除法，作(n k)(n k+1)次乘法，故消元过程中乘除法运算量为：,所以Gauss 消去法的乘除法总运算量为：,阜师院数科院第三章 线性方程组直接解法,3-17,Gauss法与Cramer法则的计算量比较,Gauss 消元法的乘除法总运算量为：,与我们曾经介绍的Cramer法则的乘除法总运算量(n21)n!+n 相比，由下表可知：当阶数越高时，Gauss消元法所需乘除法次数比Cramer法则要少得多：,Gauss 消元法的优缺点：,但其计算过程中，要求akk(k)（称为主元素）均不为零，因而适用范围小，只适用于从1到n 1阶顺序主子式均不为零的矩阵A，计算实践还表明，Gauss消元法的数值稳定性差，当出现小主元素时，会严重影响计算结果的精度，甚至导出错误的结果。,Gauss消元法简单易行。,阜师院数科院第三章 线性方程组直接解法,3-18,2 主元素法,2.1 引入主元素的必要性 对线性方程组AX=b，若其系数行列式 det(A)0，则该方程组有唯一 解，但是这一条件 不能保证所有主元素都不等于零，只要某一主元素等于零，就不能用Gauss消元法求解该方程组，即使所有主元素不等于零，但 某一主元素的绝对值很小时，Gauss消元法也是不适用的。如下例：,例2,阜师院数科院第三章 线性方程组直接解法,3-19,例2（续1）,解：为减小误差，计算过程中保留3位有效数字。按Gauss消元法步骤：,第一次消元后得同解方程组：,第二次消元后得同解方程组：,回代得解，x3=2.02，x2=2.40，x1=5.80。容易验证，方程组（3-8）的准确解为：x1=2.60，x2=1.00，x3=2.0。,显然两种结果相差很大。,阜师院数科院第三章 线性方程组直接解法,3-20,若在解方程组前，先交换方程的次序，如将（3-8）交换一行与二行改写成如下所示：,再用Gauss消元法，顺序消元后得同解方程组：,回代得解：x3=2.00，x2=1.00，x1=2.60 与准确解相同。,例2（续2）,阜师院数科院第三章 线性方程组直接解法,3-21,例2两种解法的误差分析,在例2中，对(3-8)的方程进行顺序消元时，主元a(1)11=0.50，a(2)22=0.100都比较小，以它们作 除数就增长了舍入误差,而导致计算结果不准确。,产生上述现象的原因在于舍入误差，当|a(k)kk|很小时，进行第k次消元，要用|a(k)kk|作除数，这 样就可能增大舍入误差造成溢出停机，或者导致 错误的结果。,为了在计算过程中，抑制舍入误差的增长，应尽量避免小主元的出现。如例2的第二种解法,通过交换方程次序，选取绝对值大的元素作主元 基于这种思想而导出主元素法。,阜师院数科院第三章 线性方程组直接解法,3-22,2.2 列主元素法,为简便起见，对方 程组（3-1），用 其增 广矩阵：,表示，并直接在增广矩阵上进行运算。,列主元素法的具体步骤如下：,阜师院数科院第三章 线性方程组直接解法,3-23,列主元素法,如此经过n1步，增广矩阵（3-9）被化成上三角形，然后由回代过程求解。在上述过程中，主元是按列选取的，故称为列主元法。例2中的第二种解法就是按列主元法进行的。,阜师院数科院第三章 线性方程组直接解法,3-24,2.3 全主元素法,经过nk次消元后，得到与方程组（2-1）同解的上三角形方程组，再由回代过程求解。,阜师院数科院第三章 线性方程组直接解法,3-25,主元素法举例,例6,计算过程保留三位小数。,阜师院数科院第三章 线性方程组直接解法,3-26,2.4 解三对角方程组的追赶法,在很多问题中，需要解如下形式的三对角方程组：,三对角方程组的系数矩阵为三对角阵，对于这种特殊而又简单的方程组，用前面介绍的方法求解由于有大量的零元素既占内存又浪费计算时间，显然很不经济。充分注意到三对角方程组的特点，根据顺序消元的思想导出一个简便的算法追赶法。,阜师院数科院第三章 线性方程组直接解法,3-27,追赶法的解题步骤,首先进行顺序消元，且每步将主元系数化为1，将方程组化为：,其中系数按下式计算：,然后回代求解，得：,上述追赶法能进行到底。,阜师院数科院第三章 线性方程组直接解法,3-28,追赶法举例,用追赶法解下列三对角方程组：,例4,解：首先将方程组化为（先追）：,然后回代（赶）求解：x5=0，x4=30/7，x3=6/7，x2=12/7，x1=0,可以看出，追赶法本质上还是顺序消元法，但由于计算过程中只涉及系数矩阵的非零元，因此大大节约了计算机内存与计算量，按乘除法次数进行比较，Gauss消元法约为n3/3，全主元法为n3/2，而追赶法仅为5n-3次，可见追赶法是求解三对角方程组的非常好的方法。,阜师院数科院第三章 线性方程组直接解法,3-29,3 矩阵分解法,如果用矩阵形式表示，Gauss消元法的消元过程是对方程组（3-1）的增广矩阵（A、b）进行一系列的初等行变换，将系数矩阵A化成上三角矩阵的过程，也等价于用一串初等变换阵去左乘增广矩阵，因此，消元过程可以通过矩阵运算来实现。,紧接下屏：,3.1 Gauss消元法的矩阵形式,事实上，Gauss消元法的 第一次消元相当于用初等矩阵：,阜师院数科院第三章 线性方程组直接解法,3-30,第二次消元相当于用初等矩阵：,第k次消元相当于用初等矩阵：,Gauss消元法的矩阵形式,阜师院数科院第三章 线性方程组直接解法,3-31,经过n1步消元后得到：,因为Lk(k=1,2,n1)均为非奇异阵，故它们的逆矩阵存在。容易求出：,这说明：在 的条件下，消元过程实际上是把系数矩阵A分解成单位下三角阵与上三角矩阵的乘积的过程。,Gauss消元法的矩阵形式（续）,阜师院数科院第三章 线性方程组直接解法,3-32,杜利特尔（Doolittle）分解LU分解,事实上，只要A非奇异，由上述结论，它一定可以分解成两个三角形矩阵的乘积，即：A=LU。,上述分解称为杜利特尔（Doolittle）分解，也称为LU分解，当系数矩阵完成三角分解后，对于求解方程组：Ax=b。,消元过程相当于分解A=LU及求解三角形方程组Ly=b，回代过程则是求解另一个三角形方程组Ux=y，因此，解线性方程组问题可转化为矩阵的三角分解问题。,其中：L为单位下三角矩阵，,U为上三角矩阵：,阜师院数科院第三章 线性方程组直接解法,3-33,3.2 矩阵的三角分解,正如Gauss消元法要在一定条件下才能进行到底一样，矩阵A也必须满足一定条件才能进行三角分解。,设A为n阶方阵，若A的顺序主子式Ai(i=1,2,n 1)均不为零,则矩阵A存在唯一的Doolittle分解。,定理3.1,存在性证明,唯一性证明,下面讨论如何对A进行LU分解：（1行1列）,由于两个矩阵相等就是它们的对应元素都相等，因此通过比较A与LU的对应元素，即可得到直接计算L、U的元素的公式：A=LU 即：,紧接下屏：,阜师院数科院第三章 线性方程组直接解法,3-34,对A进行LU分解,由矩阵乘法规则及比较（3-11）两端的元素，得：,即可由（3-12）求出U的第一行，L的第一列。,下面讨论一般情况，即：U的第r 行，L的第r 列：,阜师院数科院第三章 线性方程组直接解法,3-35,对A进行LU分解（i 行j 列）,一般情况下，可由aij得到计算 lij和 uij 的公式：,阜师院数科院第三章 线性方程组直接解法,3-36,计算过程应按U第1行、L第1列（第1框），U第2行、L第2列（第2框），的顺序行。,对A进行LU分解（i 行j 列）,具体分解步骤见下屏：,阜师院数科院第三章 线性方程组直接解法,3-37,对A进行LU分解的具体步骤,1.计算U的第1行，L的第1列，亦称为计算 第1框；,2.计算U的第r行，L的第r列（r=2,n），即第r框：,阜师院数科院第三章 线性方程组直接解法,3-38,矩阵A的LU分解举例,解：按分解公式（3-13），一框一框分解，每框计算时先行后列：,所以：,例5,阜师院数科院第三章 线性方程组直接解法,3-39,三角分解的紧凑格式,根据式（3-13）的特点，矩阵的三角分解可按以下格式及顺序进行。这种格式既便于记忆，又便于计算，称为紧凑格式。见下表3-1：,阜师院数科院第三章 线性方程组直接解法,3-40,（1）计算顺序：将aij，uij，lij 按表3-1列好，计算 时按框从外到内进行，每一框中先算行。从 左向右依次计算 uij；再算列，自上而下求 lij；,三角分解的紧凑格式,（2）计算方法：按行计算时，需将所求元 uij 的 对应元aij 逐次减去 uij 所在行左面各框的元 素 lij 乘以 uij 所在列上面各框相应的元 uij。按列计算 lij 时，在作上述运算后还需除以 lij 所在框的对角元 ujj。,说明：,阜师院数科院第三章 线性方程组直接解法,3-41,三角分解的紧凑格式举例,例4中矩阵A的三角分解按紧凑格式计算，结果见下表3-2：,由表可得：,阜师院数科院第三章 线性方程组直接解法,3-42,3.3 直接三角分解法,若线性方程组Ax=b的系数矩阵A完成三角分解，A=LU，那么解方程组Ax=b等价于求解两个三角形方程组Ly=b，Ux=y，即由：,再由：,可解得：,阜师院数科院第三章 线性方程组直接解法,3-43,直接三角分解法（续）,容易看出，式（3-14）与式（3-13）的运算规律相同，或者由：,可知，若将b作A的最后一列处理，在将A化为上三角阵的同时，也将b化作了y，故在利用三角分解求解方程组Ax=b时，只需将右端向量b列在表3-1的最后一列，按uij的计算方法即可求出yk（或在式（3-13）中求uij时，增加一列，j算到n+1）下表3-3是求解线性方程Ax=b的紧凑格式，其计算顺序与计算方法与三角分解相同。,按表3-3计算后，再按式（3-15）即可求出方程组的解。,阜师院数科院第三章 线性方程组直接解法,3-44,紧凑格式解线性方程组举例,解：按表3-3计算：,表3-3,例6,阜师院数科院第三章 线性方程组直接解法,3-45,所以：,紧凑格式解线性方程组举例（续）,阜师院数科院第三章 线性方程组直接解法,3-46,三角分解法的几点说明,1、用三角分解法求线性方程的乘除法运算量也是n3/3 数量级。由于在求出uij,lij和yi后，aij和bi无需保留，故上机计算时，可把L，U和y存在A，b所占的单元，回代时x取代y，整个计算过程中不需要增加新的存 贮单元。,3、完成A=LU分解后可以较容易地求出行列式|A|的值：,2、从三角分解法的推导及例中可以看出，系数矩阵的 三 角分解与右端项无关。因 而在计算多个系数矩阵 为A而右端不同的线性方程组系时，用三角分解法更 为 简便（如可用于求逆矩阵）。,阜师院数科院第三章 线性方程组直接解法,3-47,三角分解法的几点说明（续）,6、分解法的优点除上述2、3外，还有：a.可求解A2z=b,因为算A2计算量大，可用,b.可根据A的形状设计算法，当A为大型稀疏，且非零元素有规律如带状，三对角等，作分解时能充分利用A的特点，L，U能保持A的原状，即L与A的下三角相同，U与A的上三角形状相同。,4、三角分解法一般也可采用选主元的技术，以使算法更 具数值稳定性。,5、也可以把矩阵A分解成一个下三角矩阵与一个单位上三 角矩阵的乘积，矩阵的 这种分解称为克劳特（Crout）分解。,阜师院数科院第三章 线性方程组直接解法,3-48,4 平方根法与改进的平方根法,对称正定矩阵概念：,工程实际问题的计算中，线性方程组的系数矩阵常常具有对称正定性，即其各阶顺序主子式及全部特征值均大于零。矩阵的这一特性使它的三角分解也有更简单的形式，从而导出一些特殊的解法。,5.A的所有顺序主子式均为正数，即det(Ak)0,k=1,2,n)。,4.A的所有特征值0；,3.A的对角线元素aii0(i=1,2,n)；,2.A是非奇异阵，且A1也是对称正定阵；,1.A的所有顺序主子矩阵Ak(k=1,2,n)也是对称正定阵；,若A为对称正定矩阵，则：,阜师院数科院第三章 线性方程组直接解法,3-49,4.1 平方根法,定理 3.2,证明：首先A可直接作LU分解，记为A=LU1，,其中：,阜师院数科院第三章 线性方程组直接解法,3-50,定理 3.2（续）,其中U0是单位上三角阵，则由A的对称性可得：,再由唯一性得U0=LT，所以A=LDLT，并且这种分解是唯一的，又由于U1的对角线元素Ukk就是Gauss消元法的主元素：,由于LD1/2是下三角阵，因此有推论：,阜师院数科院第三章 线性方程组直接解法,3-51,Choleskg分解1,推论,若A是对称正定矩阵，则存在唯一的主对角线元素都为正的下三角阵L，使得：A=LLT,矩阵的这种分解称为Choleskg分解。用比较法可以导出L的计算公式。设：,比较A与LLT的相应元素，即由A=LLT可得：,计算顺序按列进行。,因此：,阜师院数科院第三章 线性方程组直接解法,3-52,Choleskg分解2,当完成矩阵A的Cholesky分解后，求解方 程组AX=b就转化求：,求解对称正定线性方程组的方法称为平方根法，也称为Cholesky分解法，这种方法无需选主元，计算过程也是稳定的。由于A的对称性，平方根法的乘除运算量为n3/6数量级，约为直接分解法的一半。上机计算时，所需存贮单元也少，只要存贮A的下三角部分和右端项b，计算中L存放于A的存贮单元，y,x存放在b的存贮单元。平方根法的不足之处在于需作n次开方运算。,它们的解分别为：,阜师院数科院第三章 线性方程组直接解法,3-53,平方根法举例,用平方根法解对称正定方程组：,例7,解：先分解系数矩阵A，只对A的下三角部分运算即可，所以只存放A的下三角部分：,阜师院数科院第三章 线性方程组直接解法,3-54,4.2 改进的平方根法,由定理3.2，对称正定矩阵A又可以作如下分解：,其中，L是单位下三角阵，D是对角阵，U=DLT。,由于U=DLT，即：,比较两边的对应元素可得：,首先由紧凑格式的三角分解可得：,阜师院数科院第三章 线性方程组直接解法,3-55,改进的平方根法说明,基于上述LDLT分解的方法称为改进的平方根法，其乘除运算量与Cholesky分解相当，且避免了开方运算。计算顺序按先行后列逐层分解计算；对称正定阵A完成A=ADLT=LU分解后，求解方程组：Ax=b，就转化为求解：,并且由于求y和求U的最后一列的算法完全相同，所以可将y视为U的最后一列进行计算。,按上述算法，当A为对称正定阵时，单位下三角阵L的元素不必按紧凑格式方法去求解，而只需在求得U的第k行元素后，以它们去除以ukk即得相应的L的第k列元素，这样就大大减少了计算量。,阜师院数科院第三章 线性方程组直接解法,3-56,改进的平方根法举例,用改进的平方根法求解方程组：,例8,解：对增广矩阵（A b）逐层分解可得：,则Ux=y为：,改进平方根法由于对矩阵A无正定要求（只要ukk 0），(k=1,2,n)即可，而正定要求ukk 0(k=1,2,n)，因此是求解对称方程组常用的方法。,阜师院数科院第三章 线性方程组直接解法,3-57,5 Gauss-Jordan消元法求矩阵的逆,Gauss消元法有许多变形，列主元素法是其中之一，在列主元法的基础上还可对算法进行如下的修改：在消元过程中选主元后，先将主元化为1，然后将主元所在列上，下方各元素均化为0，这样消元的结果使系数矩阵化为了单位阵，无需回代就得到了原方程之解，这种元回代过程的列主元素法称为Gauss-Jordan消元法。,Gauss-Jordan消元法比顺序消去法计算量大一点，实践中使用不多，但用它求逆阵却十分方便。因为消元过程实质上就是对系数矩阵A实行初等变换，将A化为单位阵，相当于对A左乘了一系列的初等变换阵M1，M2，Mn-1，Mn，使：,紧接下屏,阜师院数科院第三章 线性方程组直接解法,3-58,Gauss-Jordan消元法求矩阵的逆（续1）,这表明同样的一组初等变换在将A化为I 的同时，可将I 化为A1，即有：,因此，以Gauss-Jordan消元法求A的逆阵，就是要找到Mi(i=1,2,n)，以它们逐个左乘(A,I),逐列将A的对角线上的元素化为1，而其余元素化为0，最终将A化为单位阵，则I化为A的逆阵A1。,阜师院数科院第三章 线性方程组直接解法,3-59,Gauss-Jordan消元法求矩阵的逆（续2）,设增广阵为：,阜师院数科院第三章 线性方程组直接解法,3-60,这里aij(1)=a1j，上述aij(2)的计算与Gauss消元法基本上相同，仅仅由于m11与Gauss消元法中的乘数l11不相同引起第一行元素a1j(2)与aij(2)计算不相同，假若把增广阵中I的各列视为A的第n+1列，第n+2列，那么上述计算公式中的第二个下标可扩充到2n。,Gauss-Jordan消元法求矩阵的逆（续3）,阜师院数科院第三章 线性方程组直接解法,3-61,Gauss-Jordan消元法求逆阵(续4),阜师院数科院第三章 线性方程组直接解法,3-62,Gauss-Jordan消元法求逆阵(续5),阜师院数科院第三章 线性方程组直接解法,3-63,设经过k 1步后得到：,Gauss-Jordan消元法求逆阵(续6),阜师院数科院第三章 线性方程组直接解法,3-64,Gauss-Jordan消元法求逆阵(续7),其中：,阜师院数科院第三章 线性方程组直接解法,3-65,Gauss-Jordan消元法求逆阵(续8),完成n步消元后，A1放在原A的位置。,阜师院数科院第三章 线性方程组直接解法,3-66,Gauss-Jordan法求逆阵的具体步骤,按上述紧缩存贮原则，可节省存贮单元，同时还使得整个计算更简单了。可总结求逆步骤如下：,上述1，2是求第k列元素，构成Mk（即求主列）,阜师院数科院第三章 线性方程组直接解法,3-67,（计算其他元素，但少k列，k行）,用上述Gauss-Jordan法求逆阵，计算量约为n3，是Gauss消元法的3倍，为保证方法稳定性，还可选列主元，若仍按上述紧缩存贮原则，则最后需按行交换的相反次序作列交换才能得到A1。,Gauss-Jordan法求逆阵的具体步骤（续）,阜师院数科院第三章 线性方程组直接解法,3-68,Gauss-Jordan法求逆阵举例,例9,解：按紧缩存贮方式，逐次计算结果与存贮如下：,第一步：k=1，在第一列中选主元，交换1，2行，得：,k=2在第二列对角元下选主元，交换2，3行由1，2先计算第2列，由3计算其他元素（除2列2行外）而由4计算剩下的第2行的元素（这里k=2的第2列第行称为主列，主行）,第三步：k=3以a33=1/6为主元，消元后得：,交换2、3列,最后：按行交换的相反次序进行列交换：先交换2，3列，再 交换1，2列得A1。,阜师院数科院第三章 线性方程组直接解法,3-69,6 方程组的性态与条件数,无论用哪种方法求解线性方程组，一般情况下都会产生误差，本节讨论线性方程组解的误差。,方程组的解为一组数，称为解向量，近似解向量与准确解向量之差称为误差向量，为了估计误差向量的大小，以及在迭代法讨论收敛性的需要，首先需引入衡量向量与矩阵大小的度量范数。,阜师院数科院第三章 线性方程组直接解法,3-70,6.1 向量与矩阵的范数,这三个性质刻画了向量长度的基本特征，并可以用其将平面向量长度的概念推广到一般n维向量，于是有如下定义：,定义1,下屏将给出范数的种类：,阜师院数科院第三章 线性方程组直接解法,3-71,常用的向量范数,容易证明它们都满足上述三条性质。可以看出，2范数是平面向量长度计算公式在形式上的推广，也是线性代数中的内积定义。此处引入多种范数来刻画向量的大小，是为了在不同情况下用不同的范数研究问题。向量范数的证明：（只对第三条）,对范数：前面2条显然，对第三条，由于对任意实数x,y，绝对值不等式：|x+y|x|+|y|成立，因而有:,分别称为向量x的2范数，1范数，无穷范数。,阜师院数科院第三章 线性方程组直接解法,3-72,对2范数,利用实数的柯西不等式：,于是，有：,常用的向量范数（续）,阜师院数科院第三章 线性方程组直接解法,3-73,Rn中范数的等价性,例如可证明如下等价性：,所以，2范数与范数是等价的。,不难证明：,亦即1范数与范数是等价的。,事实上:Rn 中任意两种范数部是等价的。,阜师院数科院第三章 线性方程组直接解法,3-74,向量的误差,有了向量范数，就可以用它来表示方程组解向量的误差，设x是方程组Ax=b的准确解向量，是近似解向量，则:,显然，范数不同，其误差值是不一样的。,分别称为 的关于P范数的绝对误差与相对误差。,阜师院数科院第三章 线性方程组直接解法,3-75,矩阵范数,定义2,对任意n阶方阵A=(aij)nn，若对应一个非负实数|A|，满足：,则称|A|为矩阵A的范数。,与向量范数定义比较，前三条性质只是向量范数定义的推广，而第四条性质则是矩阵乘法性质的要求，它使矩阵范数在数值计算中使用更方便。,阜师院数科院第三章 线性方程组直接解法,3-76,常用的矩阵范数,常用的矩阵范数有：,它们分别叫做矩阵的范数，1范数，2范数，E范数，矩阵E范数是向量2范数的推广，矩阵范数，1范数计算容易，而矩阵2范数与ATA的特征值有关，所以又称为谱范数，它的计算较困难，但因为它有一些好的性质，所以也是常用的范数。,阜师院数科院第三章 线性方程组直接解法,3-77,常用的矩阵范数（续）,可以证明，这些范数都满足定义2。,以|A|为例，前2条性质显然成立，而对：,阜师院数科院第三章 线性方程组直接解法,3-78,最大行和矩阵范数的证明,阜师院数科院第三章 线性方程组直接解法,3-79,最大行和矩阵范数的证明（续）,阜师院数科院第三章 线性方程组直接解法,3-80,范数的相容性,在误差估计中，由于矩阵与向量会同时用到，我们总希望有上面的不等式成立，,但对任意的向量范数与矩阵范数却未必如此，因而特别地把满足此不等式的范数称为相容的，可以证明，上述常用的范数是相容的，即有：,在使用范数时，应选用相容的矩阵范数与向量范数。,分别称为 的关于P范数的绝对误差与相对误差。,有了矩阵范数，就可以用它描述矩阵的误差，设 是A的近似矩阵，称为 的残差阵，则：,阜师院数科院第三章 线性方程组直接解法,3-81,求范数举例,例10,阜师院数科院第三章 线性方程组直接解法,3-82,6.2 舍入误差的影响及算法的稳定性,前面曾介绍了多种解线性方程组的方法，由于计算机字长的限制，无论哪种方法，求解的每一步几乎都可能引入新的舍入误差，这些误差随着计算的推进而向前传播。算法不同误差的累积情况不一样，舍入误差对解的影响不大的算法称为数值稳定的算法。可以证明：主元素法、约当消去法、Cholesky分解法和追赶法都是数值稳定的算法。,阜师院数科院第三章 线性方程组直接解法,3-83,可以看出，后两个方程组与第一个方程组相比，系数矩阵或右端向量仅有0.0005以下的误差，但准确解却相差很大。,6.3 方程组的性态和条件数,数值稳定的算法是否一定能求得精度比较高的解呢？回答是不一定，解的精度还与方程组本身的性态有关，下面来考察几个例：,例11,阜师院数科院第三章 线性方程组直接解法,3-84,方程组的性态和条件数（续1）,例12,若其系数，常数项改用三位有效数字的小数表示，则方程组为：,阜师院数科院第三章 线性方程组直接解法,3-85,右端项b产生0.1%的变化引起解的变化最大变化184%。,初始数据的误差（相对）0.3%=0.003，而解的相对误差却超过50%。,例13,方程组的性态和条件数（续2）,阜师院数科院第三章 线性方程组直接解法,3-86,方程组的性态讨论 病态、良态,在许多实际问题中，线性方程组的系数矩阵和右端项的元素大多为前面计算的结果，因此上述例中的微小误差是避免不了。而对上述例中的方程组，无论用多么稳定的算法求解，计算中产生的微小误差就使解面目全非，所以这些方程组的性态是很差的。当方程组Ax=b的系数矩阵与右端向量b的微小变动（小扰动）而引起解严重失真时，称此方程组为病态方程组，其系数矩阵A称为病态矩阵，否则称为良态方程组，A称为良态矩阵，为了定量刻画方程组“病态”的程度，下面对方程组Ax=b在系数矩阵A及右端项b有扰动的几种情形进行讨论。,阜师院数科院第三章 线性方程组直接解法,3-87,此不等式表明，当右端项有扰动时，解的相对误差不超过b的相对误差的 倍。,首先考察右端项b的扰动对解的影响，设b有扰动b，A为准确,记引起解x的扰动为x,即：,方程组的性态讨论 病态、良态（续）,阜师院数科院第三章 线性方程组直接解法,3-88,方程组的性态讨论（续2）,当b为精确的而A有微小扰动A时，在A充分小时也同样可推得：,紧接下屏讨论：,阜师院数科院第三章 线性方程组直接解法,3-89,方程组的性态讨论（续3）,阜师院数科院第三章 线性方程组直接解法,3-90,方程组的性态讨论续（3）,而当A,b同时有微小扰动A，b时，则可进一步导出更一般的误差估计式：,注意到：,阜师院数科院第三章 线性方程组直接解法,3-91,在b充分小时，此式右端实际上即为：,方程组的性态讨论续（4）,阜师院数科院第三章 线性方程组直接解法,3-92,矩阵的条件数,在三种情况下得到的这三个不等式反映了解的相对误差与A及b的相对误差的关系；数|A|A-1|越小，解的相对误差也越小；反之数|A|A-1|越大，解的相对误差也越大，实际上这个数反映了解对方程组原始数据的敏感程度，揭示了矩阵A和方程组本身的性态，称之为方程组或矩阵A的条件数，记作：,cond(A)越大，A的病态程度越严重。至于cond(A)多大才算病态，这是一个相对概念，没有一个严格的数量界限。,阜师院数科院第三章 线性方程组直接解法,3-93,判断病态矩阵的几点参考,求条件数要计算逆阵的范数，很不方便，如下一些现象可作为判断病态矩阵的参考。,（1）在用主元消去法时消元过程出现小主元（如例12）,（2）矩阵A中元素间数量级相差很大；,（3）A的行列式det(A)满足（行列式值相对很大）：,（4）矩阵的某些行（或列）近似相关（如例11）。,阜师院数科院第三章 线性方程组直接解法,3-94,利用条件数判断矩阵的性态举例,A的条件数很大，所以方程组是病态的。,的特例，它是典型的“病态”阵，n越大，“病态”越严重，如n=6时，cond(A)=29106，对严重“病态”的方程组，即使用主元素法求解也难以保证数值稳定性。,阜师院数科院第三章 线性方程组直接解法,3-95,第三章,结 束,阜师院数科院第三章 线性方程组直接解法,3-96,定理3.1存在性证明,当i=1，因为A1=a11 0由3.1的讨论，存在矩阵：,假定结论对j=k1成立，即：,按归纳法原理，存在初等矩阵L1,Ln1,使得：,唯一性证明,返回,阜师院数科院第三章 线性方程组直接解法,3-97,定理3.1唯一性证明,返 回,对于唯一性，假定A有两种Dolittle分解：,由于单位下（上）三角阵的逆仍是单位下三角阵，而 为上三角阵，它们不可能相等，若要相等，必定有：,