《线性代数方程组》PPT课件.ppt

第二章,线性方程组的数值解法,确定小行星轨道,以太阳为原点在轨道平面 内建立直角坐标系,取天文测量单位,在五个不同时间观察小行星,测得坐标数据：,通过计算确定椭圆方程 a1x2 + 2a2xy + a3 y2 +2a4 x + 2a5 y + 1 = 0,a1xj2 + 2a2xjyj + a3 yj2 +2a4 xj + 2a5 yj + 1 = 0,将五个点的坐标(xj, yj) (j= 1,2,3,4,5)代入二次曲线方程,得关于a1,a2,a3,a4,a5 的方程组,在科学计算中，经常需要求解含有n个未知量 的n个方程构成的线性方程组,方程组还可以用矩阵形式表示为: Ax=b,(2.1),克莱姆法则需(n+1)(n-1) n!+n次乘除法运算 (1) 输入系数矩阵A和右端向量b;,(4) 计算并输出x1 = D1 / D;,xn=Dn/D, 结束。,(3) 对k=1,2,n用b替换A的第k列数据,并计算替换后矩阵的行列式值Dk;,(2)计算系数矩阵A的行列式值D,如果D=0,则输出错误信息结束,否则进行第(3)步;,线性方程组数值解法的分类：,线性方程组数值解法的分类：,直接法 Gauss消去法及其变形 矩阵的三角分解法,迭代法 Jacobi迭代法 Gauss-Seidel迭代法 松弛迭代,2.1 高斯消元法,1三角形方程组的解法-回代法,2.顺序Gauss（高斯）消元法是一种规则化的加减消元法。,基 本 思 想,通过逐次消元计算把需要求解的线性方程组转化成上三角形方程组，也就是把线性方程组的系数矩阵转化为上三角矩阵，从而使一般线性方程组的求解转化为等价（同解）的上三角形方程组的求解。,Gauss消元法由消元和回代两个过程组成，先讨论一个具体的线性方程组的求解。,例1 用Gauss消元法解方程组,用增广矩阵进行进算,这样，对于方程组,(2.1),用增广矩阵表示，并给出Gauss消元法的具体步骤：,或者 Ax=b,顺序Gauss消元法的消元过程可表述如下：,第一步：设 a11(1) 0 ，将第一列a11(1)以下各元素消成零,乘以矩阵A（1），b（1）的第一行再加到第i 行，得到矩阵,(i2，3，n),即依次用,其中,第二步：设 a22(2) 0 ，将第二列a22(2)以下各元素消成零,（i3，4，n）,即依次用,乘以矩阵A（2），b（2）的第二行再加到第i行，得到矩阵,其中,如此继续消元下去第n1步结束后，得到矩阵,增广矩阵A（n），b（n）对应如下上三角形方程组,这是与原线性方程组（2.1）等价的方程组.,对于等价方程组,进行回代求解，可以得到：,首先写出增广矩阵,于是，采用Gauss消元法求解方程组,(2.1),然后进行消元，采用公式,最后进行回代得到方程组的解,得到相似增广矩阵,（ik+1，k+2，n）,在编程计算时，最后的增广矩阵存放的元素是：,算法.,顺序Gauss消元法可执行的前提,定理 1 给定线性方程组 ，如果n阶方阵 的所有顺序主子式都不为零，即 则按顺序Gauss消去法所形成的各主元素 均不为零，从而Gauss 消元法可顺利执行。,注：当线性方程组的系数矩阵为对称正定或严格对角占优阵时，按Gauss消元法计算是稳定的。,例2 用Gauss消元法求解方程组：,例3 用Gauss消元法求解线性方程组,a=1 2 1 -2;2 5 3 -2;-2 -2 3 5;1 3 2 3; b=4 7 -1 0; (ab) for k=1:3 d=a(k,k); c=a(k,:); c0=b(k); for i=k+1:4 l=a(i,k)/d; a(i,:)=a(i,:)-l*c; b(i)=b(i)-l*c0; end end b(4)=b(4)/a(4,4); for k=3:-1:1 b(k)=(b(k)-a(k,k+1:4)*b(k+1:4)/a(k,k); end b=b,MATLAB 程序,ans=2 -1 2 -1 b=2 -1 2 -1,3、列主元Gauss消元法,顺序Gauss消元法计算过程中的 akk(k) 称为主元素，在第k步消元时要用它作除数,则可能会出现以下几种情况,若出现 akk(k) 0，消元过程就不能进行下去。,akk(k) 0 ，消元过程能够进行，但若|akk(k)| 过小，也会造成舍入误差积累很大导致计算解的精度下降。,例4: 在四位十进制的限制下，试用顺序Gauss消元法求解如下方程组,此方程组具有四位有效数字的精确解为,x117.46，x245.76，x35.546,解 用顺序Gauss消元法求解，消元过程如下,经回代求解得 x35.546，x2100.0，x1104.0,和此方程组的精确解相比,x35.546 ，x245.76， x117.46,有较大的误差。,对于此例，由于顺序Gauss消元法中的主元素绝对值非常小，使消元乘数绝对值非常大，计算过程中出现大数吃掉小数现象，产生较大的舍入误差，最终导致计算解 x1104.0 和 x2100.0 已完全失真。,为避免这种现象发生，可以对原方程组作等价变换，再利用顺序Gauss消元法求解。,写出原方程组的增广矩阵：,针对第一列找出绝对值最大的元素，进行等价变换：,求得方程的解为：x35.546，x245.76，x117.46,精确解为： x35.546 ，x245.76， x117.46,列主元Gauss消元法与顺序Gauss消元法的不同之处在于：,后者是按自然顺序取主元素进行消元,前者在每步消元之前先选取主元素然后再进行消元,下面将列主元Gauss消元法的计算步骤叙述如下：,给定线性方程组 Axb, 记A（1）, b（1） A,b，列主元Gauss消去法的具体过程如下：,1. 首先在增广矩阵A（1）,b（1）第一列的n个元素中选取绝对值最大的一个作为主元素，并把此主元素所在的行与第一行交换，即,2. 其次进行第一步消元得到增广矩阵A（2）,b（2），在矩阵A（2）,b（2） 第二列的后 n1个元素中选取绝对值最大的一个作为主元素，并把此主元素所在的行与第二行交换，即,3. 再进行第二步消元得到增广矩阵A（3），b（3）。按此方法继续进行下去，经过n1步选主元和消元运算，得到增广矩阵A（n），b（n），它对应的方程组,A（n）xb（n）,是一个与原方程组等价的上三角形方程组，可进行回代求解。,容易证明，只要det（A）0，列主元Gauss消去法就可以顺利完成，即不会出现主元素为零或者绝对值太小的情形出现。,选列主元过程：,一、求主元 alk 使得 | alk | =max |akk|, |ak+1,k|, , |ank | ;,二、判断：若 |alk| ，则输出错误信息并停机，否则 转三；,三、判断：若 l k ，则交换增广矩阵中第 l 行与 k 行 的元素，否则不交换。,例5 用列主元法解,第一列中绝对值最大是10，取10为主元；,n 阶方程组第 k 轮消元时，选第 k 列的后 (n-k+1) 个元素中绝对值最大者作主元。,x1 = 0 x2 = -1 x3 = 1,第二列的后两个数中选出主元 2.5；,高斯消元法是一种顺序消元法。消元过程按方程和未知数的顺序依次作消元计算，缺乏一定的灵活性。计算中有除数为零时，算法便无法实现；即使所有的除数都不为零，可以计算出方程组的解，也不能保证计算结果有很高的精度。,m21 = 2/0.00001 = 2 105,= 0.000003 106 0.40000 106 = -0.40000 106,= 0.000002 106 0.20000 106 = -0.20000 106,x2 = 0.50000, x1 = 0.00000.,解法一：高斯消元法,m21 = 0.00001/2 = 0.5 10-5,=2 - (0.5 10-5) 3,x2 = 0.50000, x1 = 0.25000.,解法二：,列主元消元法,= 0.20000 10 0.0000015 10 = 0.20000 10,=1 - (0.5 10-5) 2,= 0.10000 10 0.000001 10 = 0.10000 10,2.2 矩阵的直接分解法, 高斯消元法的矩阵形式,每一步消去过程相当于左乘初等下三角矩阵Lk,A 的 LU 分解 ( LU factorization ),第一个方程组的系数矩阵为下三角矩阵，第二个方程组的系数矩阵为上三角矩阵，两个方程组都非常容易求解。,将A=LU 称为矩阵A的三角分解,这时线性方程组为：,令,则有,定理2：（矩阵的三角分解）设A为n n实矩阵，如果 解AX = b用高斯消去法能够完成（限制不进行行的交 换，即 ），则矩阵A可分解 为单位下三角矩阵L与上三角知阵U的乘积。 A = LU 且这种分解是唯一的。,注： （1） L 为单位下三角阵而 U 为一般上三角阵的分解称为多利特尔（Doolittle） 分解； （2）L 为一般下三角阵而 U 为单位上三角阵的分解称为库朗（Courant ）分解。, Doolittle分解法 ：,根据 A=LU 有等式成立：,比较等式两端对应元素，有,可以解得：,当 i=1 时,当 j=1 时,当 i 1 时,当 j 1 时,下面，我们对具体矩阵进行Doolittle 三角分解。,为了表示和存储方便，可以将分解后的两个矩阵用一个矩阵表示,例7 利用Doolittle三角分解法分解矩阵,解：分解时用到如下公式,1,2,3,4,1,1,1,2,6,12,3,7,6,24,6,24,1,2,3,4,1,1,1,2,6,12,3,7,6,24,6,24,可以写成：,这时，矩阵的三角分解,如果我们要求解方程组,则由,得到,例8：利用Doolittle三角分解方法解线性方程组,解: 进行三角分解ALU，可以对增广矩阵A,b作 三角分解:,1 2 3 -2,-3,2,2,-3,-1,3,3,17,得到,1 2 3 -2,-3,2,2,-3,-1,3,3,17,这时，相应的方程组为：,x135,x28，,例9： 用矩阵分解方法解 AX = b,LY = b Ux = Y,y1 = 6 y2 = -4 y3 = -4,x1 = -13 x2 = 8 x3 = 2,A=LU,不选主元的LU分解,A=2，3，4；3，5，2；4，3，30 n=max(size(A); for k=1:n-1 d=A(k,k); for i=k+1:n m=A(i,k)/d; for j=k+1:n A(i,j)=A(i,j)-m*A(k,j); end A(i,k)=m; end end,A = 2.00 3.00 4.00 1.50 0.50 -4.00 2.00 -6.00 -2.00, 平方根法 ：,在实际应用中，常见一类非常重要的线性方程组Axb，其中A为对称正定矩阵，即A是对称的且对任何非零向量 x 都有xTAx0。本节将对这类方程组导出更有效地三角分解求解方法，称之为平方根法。,设A为对称正定矩阵，那么A的所有顺序主子式均大于零，根据定理2.2，存在惟一三角分解 ALU，即,记 Ak（1kn）为A的 k 阶顺序主子阵，则det（Ak）为A的k阶顺序主子式。由上式，利用矩阵分块运算规则，容易验证,det（Ak）u11 u22 ukk,那么由 det（Ak）0，可知 ukk0，k1，2，n,这时，将上面的矩阵表示为：,即： A=LDM ，其中 DM=U， M=D-1U。,当AAT为对称矩阵时，根据 ALDM, 得到,ATMT D LT,再根据矩阵三角分解的唯一性,可知 MLT。于是,ALDLT,则有,令,如果对称正定矩阵A具有如下分解 AGGT，其中G为下三角矩阵，则称其为对称正定矩阵的 Cholesky（乔列斯基）分解。,为表示方便，可以记,给定对称正定方程组 Axb，对 A 进行 Cholesky分解ALLT，则原方程组等价于 LLTxb,Lyb LTx y,解此方程组即可得到原方程组的解x ，这就是求解方程组的平方根法。,下面，我们通过比较矩阵的对应元素给出对称正定矩阵的平方根分解法。已知,即： LLTxb，等价于,比较对应元素：,当 i = j 时,当 j i 时,解得,于是，根据计算公式,可以对对称正定矩阵进行平方根分解,l11,l21,l22,l31,l32,l33,ln1,ln2,ln3,lnn,关于方程组 Ax=b , 如果对系数矩阵进行了平方根分解 ALLT，则将方程组化为： Lyb , LTxy,解得,于是，关于系数矩阵是对称正定矩阵的线性方程组Ax=b的求解，分两步进行：,第一步：系数矩阵的平方根分解,第二步：解等价方程组,例10 用平方根法求解对称正定方程组,解: 首先进行A 的Cholesky 分解,ALLT,2,-0.5,0.5,2,1.5,1,得 y12，y23.5 ，y31,得 x11，x21，x31,求解Lyb:,再求解 LTxy:, 解三对角方程组的追赶法,考虑如下形式的三对角方程组问题,其中，系数矩阵满足条件:,该三对角方程组的增广矩阵 为:,第一步：第一行元素除以b1，取 1 = c1 / b1，y1 = f1 / b1，可得 增广矩阵 ：,第二步：从出发，用初等变换作 n-1 轮消元。作第 k 轮消元时，将矩阵中第 k 行元素乘以 - ak+1 加到第 k+1 行元素上，然后将第 k+1 行主元单位化( k = 1 , 2, , n-1 )。最后的增广矩阵 ：,在消元过程中所用除法次数为 2n - 1，所用乘法次数为 2(n - 1)。,根据矩阵初等变换的性质知，原三对角方程组,等价于如下方程组,对这一特殊的上三角方程组，在回代过程中只需用到 n 1 次乘法就可以求出方程组的解。,追赶法算法框图,例11 用追赶法求解方程组,最后得增广矩阵,所以，方程组的解为,a=0 1 1 1 1; b=4 4 4 4 4; c=1 1 1 1 0; f=2 1 1 1 2; d=b(1); f(1)=f(1)/d; for k=2:5, c(k-1)=c(k-1)/d; d=b(k)-a(k)*c(k-1); f(k)=(f(k)-a(k)*f(k-1)/d; end clear d a b for k=4:-1:1, f(k)=f(k)-c(k)*f(k+1); end f,f = 0.4808 0.0769 0.2115 0.0769 0.4808,MATLAB 程序,