1第五章线性代数方程组的数值解法线性方程求解问题是科学研究和工程计算中最常见的问题。如电学中的网络问题、工程力学中求解连续力学体(微分方程)问题的差分方法、有限元法、边界元法及函数的样条插值、最小二乘拟合等,都包含了解线性方程组问题。因此,线性方程组的解法在数值计算中占有极其重要的地位。对于n阶线性方程组Axb,若det()0A,则方程组有惟一解。由克莱姆(Cramer)法则,其解为det()(1,2,,)det()iiAxinA,其中iA为用向量b代替A中第i列向量所得矩阵。每个n阶行列式共有!n项,每项都有n个因子,所以计算一个n阶行列式需做(1)!nn次乘法,我们共需要计算1n个行列式,要计算出ix,还要做n次除法,因此用Cramer法则求解要做2(1)!nnn次乘除法(不计加减法),计算量十分惊人。如30n时,就需作约352.3810次乘法。可见Cramer法则在理论上是绝对正确的,但当n较大时,在2实际计算中却是不可行的。因此寻求有效的数值计算方法就成为非常必要的课题。线性方程组的类型很多,若按其系数矩阵阶数的高低和含零元素多少,大致可分为两类:一类是低阶稠密线性方程组,即系数矩阵阶数不高,含零元素很少。另一类是高阶稀疏线性方程组,即系数矩阵阶数高,零元素占绝对优势(比如占70%以上)。线性方程组的数值解法也可分为两大类:直接法和迭代法。直接法是在没有舍入误差的情况下,通过有限步运算可以得到方程组精确解的方法。但是,在实际计算时,由于初始数据变为机器数而产生的误差以及计算过程中所产生的舍入误差等都要对解的精确度产生影响,因此直接法实际上也只能算出方程真解的近似值。常用的有效算法是Gauss消去法和矩阵的三角分解法。迭代法是用某种极限过程去逼近准确解的方法。如对任意给定的初始近似解向量(0)x,按照某种方法逐步生成近似解序列(0)(1)(),,,,,kxxx3使极限()*limkkxx为方程组的解。因为在实际计算时,只能做到有限步,所以得到的也是近似解。常用的迭代法主要有Jacobi迭代法、Gauss-Seidel迭代法、逐次超松弛法及共轭斜量法等。直接法的优点是运算次数固定,并且可以事先估计。它的缺点是通常需要存储系数矩阵和常数项的所有元素,因而所需存贮单元较多,编写程序较复杂,一般适用于求解低阶线性方程组;迭代法的优点是原始系数矩阵始终不变,且只需存储原方程组系数矩阵中的非零元素,因而算法简单,编写程序较方便,且所需存贮单元也较少,所以被广泛用于求解高阶稀疏线性方程组。缺点是存在收敛性和收敛速度问题,因而只能对具有某些性质的方程组才适用。1消元法1.1三角方程组的解法形如4111122111,111,1nnnnnnnnnnnnnuxuxuxyuxuxyuxy(1.1)的方程组称为上三角形方程组。写成矩阵形式111,11111,11,11nnnnnnnnnnnnuuuxyuuxyuxy简写为Uxy。若det()0U,则(1.1)有惟一解1(1,,1)nnnnnkkkjjkkjkxyuknxyuxu我们称求解上三角方程组(1.1)的过程为回代过程。同时也看到求解这类方程组是件容易的事,所以对一般形式的方程组,应设法将它化为(1.1)式的形式,然后再求解。1.2Gauss消去法考虑方程组5111122111122nnnnnnnnaxaxaxbaxaxaxb(1.2)其矩阵形式为Ax=b其中111111nnnnnnaaxbaaxbAxb,,化线性方程组(1.2)为等价的三角形方程组的方法有多种,由此导出不同的直接方法,其中Gauss消去法是最基本的一种方法。Gauss消去法分消元计算和回代求解两个过程。为后面的符号统一起见,记方程组(1.2)为(1)(1)Axb其中(1)(1)(1)()(),ijijaaAbb。消元计算第一步:就是要将(1)A的第一列主对角元以下的元素全约化为0。设(1)110a,计算(1)(1)1111(2,3,,)iimaain6用1im乘(1.2)的第1个方程,加到第i个方程上,完成第一步消元,得(1.2)的同解方程组(1)(1)(1)(1)1111211(2)(2)(2)22222(2)(2)(2)200nnnnnnnxaaabxaabxaab(1.3)简记为(2)(2)Axb其中(2)(2),Ab的元素的计算公式为(2)(1)(1)11(2)(1)(1)11,,2,,,2,,ijijijiiiaamaijnbbmbin假设前1k步消元完成后,得(1.2)的同解方程组为(1)(1)(1)(1)(1)11112111(2)(2)(2)(2)222222()()()()()()knknkkkkkkknkkkknnknnnxaaaabxaaabxaabxaab(1.4)简记为()()kkAxb。第k步:就是要将()kA的第k列主对角元以下7的元素全约化为0。设()0kkka,计算()()(1,,)kkikikkkmaaikn用ikm乘(1.4)的第k个方程加到第i个(1,ik,)n方程上,完成第k步消元。得同解方程组(1)(1)kkAxb其中(1)(1),kkAb元素的计算公式为(1)()()(1)()()(,1,,)kkkijijikkjkkkiiikkaamaijknbbmb继续上述过程,完成1n步消元后,(1.2)化成同解的上三角方程组(1)(1)(1)(1)1111211(2)(2)(2)22222()()nnnnnnnnxaaabxaabxab(1.5)简记为()()nnAxb。回代求解因det()0A,故()0nnna,于是()()()()()1(1,,1)nnnnnnnkkkkkkllkklkxbaxbaxakn(1.6)8Gauss消去步骤能顺利进行的条件是()0kkka,1,,1kn,现在的问题是矩阵A应具有什么性质,才能保证此条件成立。若用iD表示A的顺序主子式,即11111,,iiiiiaaDinaa,有下面定理定理1约化的主元素()0(1,,)iiiaik的充要条件是矩阵A的顺序主子式01,,iDik,。证明必要性。因主元素()0(1,,)iiiaik,可进行1()mmk步消元,每步消元过程不改变顺序主子式的值,于是(1)()11,,0()mmmmDaamk必要性得证。用归纳法证明充分性。1k时命题显然成立。假设命题对1k成立。设01,,iDik,。由归纳法假设有()0(1,,1)iiiaik,Gauss消去法可以进行1k步,A约化为9()()()1112()220kkkkAAAA其中()11kA是对角元为(1)(2)(1)11221,1,,,kkkaaa的上三角阵。因()kA是通过1k步消去法得到的,每步消元过程不改变顺序主子式的值,所以A的k阶顺序主子式等于()kA的,即()(1)(1)()11111,1()*det,,0kkkkkkkkkkkADaaaa由0kD知(),0kkka,充分性得证。1.3Gauss消去法的计算量1)消元过程的工作量第(1,,1)kkn步消元过程:计算乘子()(),1,,kkikikkkmaaikn需作nk次除法运算;第k行乘ikm加到第i行,需要乘法和减法各1nk次(第k列元素无需计算),共nk行(1,,ikn),所以需要乘法()(1)nknk次,故消元过程中的乘、除和减法运算量为乘除法次数:10121[()()(1)](1)(1)23nknnnknknknn加减法次数:121()(1)(1)3nknnknkn2)回代过程的工作量计算()()()1(1,,1)nkkkkkkllkklkxbaxakn需要(1)nk次乘除法和()nk次加减法运算,整个回代过程的运算量为:乘除法:11(1)(1)2nknnkn减法:11()(1)2nknnkn所以Gauss消去法的总运算量为乘除法:3322(1)(1)(1)322333nnnnnnnnnn(当n较大时)加减法:32325(1)(1)323263nnnnnnnn(当n较大时)当10n时,用Gauss消去法约需430次乘除运算,而用Gramer法则约需210!10次乘除法。111.4Gauss消去法的矩阵形式如果用矩阵形式表示,Gauss消去法的消元过程是对方程组(1.2)的增广矩阵[]A,b进行一系列初等变换,将系数矩阵A化成上三角形矩阵的过程,也等价于用一串初等矩阵去左乘增广矩阵,因此消元过程可以通过矩阵运算来表示。第一次消元等价于用初等矩阵21111111nmLm左乘矩阵(1)(1)[,]Ab,即(2)(2)(1)(1)1[,][,]LAbAb一般地,第k次消元等价于用初等矩阵1,,1111kkknkLmm12左乘矩阵()()[,]kkAb,即(1)(1)()()1[,][,]kkkkLAbAb经过1n次消元后得到(1)(1)()()121[,][,]nnnnLLLAbAb有(1)()121(1)()121nnnnnnLLLAALLLbb将上三角矩阵()nA记为U,得到111121nLLLAULU其中211113212112,11111nnnnnmmLLLmmmL为单位下三角阵。这说明,消元过程实际上是把系数矩阵A分解为单位下三角阵与上三角矩阵的乘积的过程。这种分解称为杜利特尔(Doolittle)分解,也称为LU分解。定理2(矩阵的LU分解)设A为n阶方阵,13如果A的顺序主子式0(1,,1)iDin,则存在惟一的单位下三角阵L和上三角阵U,使ALU。证明以上的分析已证明了A可作LU分解。下面证明分解的惟一性。如果A非奇异,设A有两个分解式11ALULU其中1,LL都是单位下三角阵,1,UU都是上三角阵。因A非奇异,所以11,,,LLUU都可逆。于是1111LLUU左边为单位下三角阵,右边为上三角阵,所以1111LLUUI即有11,LLUU,惟一性得证。若A为奇异阵,因其1n阶顺序主子式不等于零,故可记12112330,,10nnnnAALUUALUAaLu其中1A为1n阶非奇异阵,1L为1n阶单位下三角阵,由ALU,得11121233132,,,nnnnALUALUALUaLUu14由已证明的结果知111ALU的分解是惟一。所以11212331,ULALAU亦是惟一确定的。进而,nnu也是惟一确定的。定理2的逆命题是:定理3设n阶方阵A非奇异,若A有惟一的分解ALU(L为单位下三角阵,U为上三角阵),试证A的顺序主子式0,1,2,,1iDin。(证明留给读者)2主元素法前面已指出Gauss消去法必须在各个约化主元素()0(1,,1)kkkakn下才能进行。现在还需