解线性方程组的直接法.ppt

1,主讲教师：高小辉E-mail:fzlcstar,第四章线性方程组AX=B的数值解法(TheSolutionofLinearSystemsAX=B),2,求解,4.1引言,许多实际问题可归结为线性（代数）方程组,机械设备、土建结构的受力分析；经济计划,输电网络、管道系统的参数计算；企业管理,大型的方程组需要有效的数值解法。,数值解法的稳定性和收敛性问题需要注意。,3,小行星轨道计算问题,4,天文学家要确定一小行星的轨道,在轨道平面建立以太阳为原点的直角坐标系.在坐标轴上取天文单位(地球到太阳的平均距离),对小行星作5次观察,测得坐标数据(x1,y1),(x2,y2),(x3,y3),(x4,y4),(x5,y5),将数据代入椭圆的一般方程:a1x2+a2xy+a3y2+a4x+a5y+1=0,得a1xk2+a2xkyk+a3yk2+a4xk+a5yk+1=0(k=1,2,3,4,5),5,即求解方程组,可确定椭圆方程(小行星轨道方程),6,对一般线性方程组:AX=b,其中,当系数矩阵A的行列式|A|0时,则方程组有唯一解.,7,求解线性方程组:AX=b的一般过程:,输入：A，b,解方程组算法,输出：X,直接法：经过有限步算术运算求得精确解迭代法：从初始解出发，逐步求出近似解来逼近,在求解小型(未知数较少)方程组时,直接法很有效.在求解大型方程组时,迭代法是最有效的方法.,8,4.2高斯消去法(GaussianElimination),1.上三角线性方程组(Upper-triangularLinearSystem)定义4.1NN矩阵A=aij中的元素满足对所有ij,有aij=0,则称NN矩阵A=aij为上三角矩阵。如果A中的元素满足对所有ij,有aij=0,则称NN矩阵A=aij为下三角矩阵。,4.2.1顺序高斯消去法,9,10,AX=B上三角线性方程组表示为：a11x1+a12x2+a13x3+a1n-1xn-1+a1nxn=b1a22x2+a23x3+a2n-1xn-1+a2nxn=b2a33x3+a3n-1xn-1+a3nxn=b3.an-1n-1xn-1+an-1nxn=bn-1annxn=bn,2.回代(BackSubstitution),设AX=B是上三角线性方程组，如果：akk0,k=1,2.n,则方程组存在唯一解。,11,证明：,a11x1+a12x2+a13x3+a1n-1xn-1+a1nxn=b1a22x2+a23x3+a2n-1xn-1+a2nxn=b2a33x3+a3n-1xn-1+a3nxn=b3.an-1n-1xn-1+an-1nxn=bn-1annxn=bn,唯一,用归纳法可证明xn-1,xn-2.x1是唯一的,12,例4.2：证明下列线性方程组无解4x1x2+2x3+3x4=207x3-4x4=-76x3+5x4=43x4=6,例4.1：利用回代法求解线性方程组4x1x2+2x3+3x4=202x2+7x3-4x4=-76x3+5x4=43x4=6,x4=6/3=2,x3=(4-5x4)/6=-1,x2=(-7-7x3+4x4)/-2=-4,x1=(20+x2-2x3-3x4)/4=3,x4=2,x3=-1,x3=1/7,a22=0,0 x2+,13,例4.3：证明下列线性方程组有无穷解4x1x2+2x3+3x4=200 x2+7x3-0 x4=-76x3+5x4=43x4=6,x4=2,x3=-1,x3=-1,x2=4x1-16,无穷解,a22=0,14,复习,如A是NN矩阵,则A的行列式为：,划掉A的第i行和第j列后的行列式,第i行展开,第j列展开,15,16,3.如果NN矩阵A=aij是上三角矩阵或下三角矩阵，则：,定理：A是NN方阵，线性方程组AX=B有唯一解当且仅当det(A)0,17,高斯消元法：,18,4初等变换(ElementaryTransformation),下列三种变换可使一个线性方程组变换成另一个等价的线性方程组交换变换：对调方程组的两行比例变换：用非零常数乘方程组的某一行替换变换：将方程组的某一行乘一个常数再加到另一行,19,例4.4：求抛物线y=A+Bx+Cx2,它经过三点(1,1),(2,-1),(3,1)。,(2)-(1),(3)-2(2),C=2,B=-8,A=7,抛物线方程y=7-8x+2x2,A+B+C=1(1)A+2B+4C=-1(2)A+3B+9C=1(3),(3)-(1),(3-1),(3-2),20,上述求解方程组的方法就是高斯(Gauss)消元法。从式(4-1)到(4-2)的过程称为消元过程而由(4-2)求出C、B、A的过程称为回代过程。因此用高斯消去法求解性方程组要经过消元和回代两个过程。,21,可将线性方程组AX=B的系数存在N(N+1)增广矩阵中,22,5初等行变换(ElementaryRowOperation),对增广矩阵进行如下变换可得到一个等价的线性方程组交换行变换：对调矩阵两行比例行变换：用非零常数乘矩阵某一行所有的元素替换行变换：将矩阵的某一行的所有元素乘一个常数再加到另一行对应的元素上去。,23,4.2.2主元素（Pivot）高斯消去法系数矩阵A中的元素arr用来消去akr(k=r+1,r+2,N),arr为第r个主元，第r行为主元行。,例4.6用增广矩阵表示下列线性方程组，并求等价上三角线性方程组和方程组的解。,x1+2x2+x3+4x4=132x1+0 x2+4x3+3x4=284x1+2x2+2x3+x4=20-3x1+x2+3x3+2x4=6,24,x1+2x2+x3+4x4=132x1+0 x2+4x3+3x4=284x1+2x2+2x3+x4=20-3x1+x2+3x3+2x4=6,0761445,-主元行2,-主元行4,-主元行-3,121413,042-52,06215-32,25,-主元行-1.75,-主元行1.5,009.55.2548.5,121413,042-52,0057.5-35,00009-18,-主元行1.9,回代：x4=2,x3=4,x2=-1,x1=3,消元过程,26,有回代的高斯消去法(GaussianEliminationwithBackSubstitution),如果A是NN非奇异矩阵（存在A-1），则存在线性方程组UX=Y与线性方程组AX=B等价，这里U是上三角矩阵，并且akk0。当构造出U和Y后，可用回代法求解UX=Y，并得到方程组的解X。,27,消元,记,28,Step1：设，计算因子,将增广矩阵第i行mi1第1行，得到,其中,29,一般地，假定已完成了(k-1)步消元，即已将A(1)|B(1)转化为以下形式：,30,Stepk：设，计算因子,共进行？步,n1,且计算,31,回代,32,消元法步骤(以4阶为例)：,增广矩阵,33,计算3个消元因子(乘子向量)(1)设a110,m21m31m41T=a21a31a41T/a11,用-m21乘矩阵第一行后加到矩阵第二行;用-m31乘矩阵第一行后加到矩阵第三行;用-m41乘矩阵第一行后加到矩阵第四行;,34,用-m31乘矩阵第一行后加到矩阵第三行;,35,用-m41乘矩阵第一行后加到矩阵第四行;,36,第二轮消元后,(2)设a22(1)0,计算第二个消元因子m32m42T=a32(1)a42(1)T/a22(1),用-m32乘矩阵第二行后加到矩阵第三行;用-m42乘矩阵第二行后加到矩阵第四行;,37,第三轮消元后,用-m43乘矩阵第三行后加到矩阵第四行;,(3)设a33(2)0,计算第三个消元因子:m43=a43(2)/a33(2),对应三角形方程组,38,消元过程中,消元因子可排列为一矩阵,用回代算法即可得出方程组的解.,39,由初等变换的知识可得.LU=A,在消元过程中,我们得到两个矩阵,此即为矩阵A的LU分解.,40,高斯消元法的主要缺陷1.计算消元因子mik=aik/akk当分母akk为零时,计算无法进行;,2.分母akk的绝对值非常小时,也将会对计算结果误差产生很大影响,解决办法,41,选主元以避免app(p)=0如果app(p)=0,寻找第p行下满足akp(p)0的第一行，设为第k行，然后交换行k和行p，使新app(p)0。如果app(p)0则不交换。此方法称为平凡选主元方法。,42,将列主元所在行与第k行交换后,再按上面的高斯消元法进行下去,此法即称为列主元消元法。,选|aik(k)|(i=k,n)最大的一个(列主元),选主元以减少误差,在每一步消元过程中取系数子矩阵的第一列元素中绝对值最大者作主元。,43,例4.7:值x1=x2=1.000是如下方程组的解：1.133x1+5.281x2=6.41424.14x1-1.210 x2=22.93用4位有效数字及选主元的高斯法求方程的解。行2-行124.14/1.133:1.133x1+5.281x2=6.414-113.7x2=-113.8,x2=1.001,x1=0.9956,44,例4.8：24.14x1-1.210 x2=22.931.133x1+5.281x2=6.414行2-行11.133/24.14:24.14x1-1.210 x2=22.935.338x2=5.338,x2=1.000,x1=1.000,45,例用列主元法解,第一列中绝对值最大是10，取10为主元,46,第二列的后两个数中选出主元2.5,47,X1=0X2=-1X3=1,N阶方程组第k轮消元时，选第k列的后(n-k+1)个元素中绝对值最大者作主元。,48,4.3矩阵三角分解法（TrianguarFactorization）定义如果非奇异矩阵A可表示为下三角矩阵L和上三角矩阵U的乘积A=LU，则A存在一个三角分解。4.3.1高斯消去法和矩阵三角分解法,a11x1+a12x2+a13x3+a1n-1xn-1+a1nxn=b1a21x1+a22x2+a23x3+a2n-1xn-1+a2nxn=b2a31x1+a32x2+a33x3+a3n-1xn-1+a3nxn=b3an1x1+an2x2+an3x3+ann-1xn-1+annxn=bn,49,=,50,51,UX=Y,LY=b,Y,X,52,由初等变换的知识可得.LU=A,在高斯消元过程中,我们得到两个矩阵,此即为矩阵A的LU分解.,53,例4.3.2求解：,x1+2x2+4x3+x4=212x1+8x2+6x3+4x4=523x1+10 x2+8x3+8x4=794x1+12x2+10 x3+6x4=82,54,y1=21y2=10y3=6y4=-24,55,x4=4x3=3x2=2x1=1,56,4.3.2直接三角分解法例4.3.1使用高斯消去构造下列矩阵的三角分解：,57,定理,设无行交换变换的高斯消去法可求解一般线性方程组AX=B，则矩阵A可分解为一个下三角矩阵L和一个上三角矩阵U的乘积：A=LU，而且L对角线元素为1，U的对角线元素非零。得到L和U后，可通过如下步骤得到X：1、利用前向替换法对方程组LY=B求解Y2、利用回代法对方程组UX=Y求解X,58,4.4置换矩阵可能一个非奇异矩阵不能直接分解为A=LU,例4.4.1证明下列矩阵不能直接分解为A=LU：,59,=,1=u114=m21u11=m21-2=m31u11=m31,2=u128=m21u12+u22=4*2+03=m31u12+m32u22=(-2)*2+m32*0=-4,60,定义,一个NN置换矩阵P是一个在每一行和每一列只有一个元素为1，而其它元素为0的矩阵。P的列是单位矩阵行的置换，可表示为：P=Ek1Ek2Ekn,61,定理,设P=Ek1Ek2Ekn是一个置换矩阵。PA是一个新的矩阵，它的行是将A中的行按行k1A,行k2A，行knA调整顺序后形成的。,例4.4.2设A为44矩阵，P为上式中的置换矩阵，则PA矩阵中的行是将A中的行调整顺序后形成的，顺序为第1，2，3，4行对应于A中的第2，1，4，3行。,62,=,63,定理,如果A是非奇异矩阵，则存在一个置换矩阵P，使得PA存在三角分解：PA=LU,例4.4.3设A为非奇异矩阵,=,64,4.5扩展高斯消去过程,定理,（非直接分解：PA=LU）,设A是一NN矩阵。假设高斯消去法可求解经过行变换的一般线性方程组AX=B。则存在一个置换矩阵P，使得PA可分解为一个下三角矩阵L和一个上三角矩阵U：PA=LU。而且L的主对角线元素为1，U的主对角线元素非零。可用如下4步求出X：1、构造矩阵L，U，P2、计算列向量PB3、用前向替换法对方程组LY=PB求解Y4、用回代法对方程组UX=Y求解X,如何求P,65,通过选列主元,PAX=PB,AX=B,66,x1+2x2+4x3+x4=212x1+8x2+6x3+4x4=523x1+10 x2+8x3+8x4=794x1+12x2+10 x3+6x4=82,A=,1、构造矩阵L，U，P,67,=,U,2、计算列向量PB,L,P,=,68,3、用前向替换法对方程组LY=PB求解Y,=,4、用回代法对方程组UX=Y求解X,=,69,%构造P,p=zeros(n,n);fori=1:np(i,i)=1;endfork=1:n-1amax,j=max(abs(a(k:n,k);c=a(k,:);a(k,:)=a(j+k-1,:);a(j+k-1,:)=c;c=p(k,:);p(k,:)