Gauss消元法解解线性方程组

摘要本文叙述了Gauss顺序消元法解线性方程的算法思想以及其求解过程，同时简要叙述了Gauss主元素消元法以及Gauss全主元消元法。紧接着给出了GaussSeidel迭代法的算法思想，本文给出了这三个消元方法以及一个迭代法的算法流程图，由于全主元消元法是前两个算法的基础上改进而来，故本文采用第三种方法进行编程计算，前两种方法不再重复编程，然后给出一个实例的计算结果，运行时间，在文章最后分析该实例的计算结果，针对同一实例，又采用GaussSeidel方法编程实现，然后对结果进行分析和对比。最后给出了本人在编程时遇到的一些问题和解决办法。关键词：Gauss顺序消元法Gauss主元素消元法Gauss全主元消元法一、算法的简要描述1.1Gauss顺序消元法Gauss消元法在中学里已经学习过，其方法实质，就是运用初等变换，将线性方程组Axb转化为同解的上三角矩阵方程组1UxLb（1.1.1）其中，U为上三角矩阵，L为下三角矩阵。然后对式（1.1.1）进行回代求解，即得方程组的解。手算的过程是非常清楚的，现在需回答的是计算机求解，如何实现上述计算过程。设线性方程组为1111221331121122223322112233nnnnnnnnnnnaxaxaxaxbaxaxaxaxbaxaxaxaxb写成矩阵形式为1112111212222221222mmmnnaaaxbaaaxbaaaxb（1.1.2）设线性方程组如上式所示，记(1)AA，(1)bb，与是增广矩阵具有形式(1)(1)[][]AbAb，此时方程组为(1)(1)Axb。第一次消元。设(1)110a，为将第二个方程至第n个方程的1x系数(1)1ia消成零，构造乘数(1)11(1)11iiala(2,3,,)in用1il乘以矩阵(1)[Ab]（1）的第一行加到第i行上(2,3,,)in得(1)(1)(1)(1)(1)11121311(2)(2)(2)(2)(2)222322(2)(2)(2)(2)230[Ab]0nnnnnnnaaaabaaabaaab（2）其中，(2)(1)(1)11(2)(1)(1)11aa(,2,3,,)(,2,3,,)ijijijiiilaijnbblbijn假设经过1k次消元得同解方程组()()kkAxb，此时(1)(1)(1)(1)(1)11121311(2)(2)(2)(2)k()222322()()()()()nknkkknkkkknknnnaaaabaaabAbabaab（）(1.1.3)假若()0kkka，则第k次消元如下：构造乘数()()(1,,)kikikkkkalikna用ikl乘以矩阵()()kkAb第k行加到第i行上(1,,)ikn得同解方程组(1)(1)kkAxb，其中，(1)(1)[]kkAb中的元素计算公式为(1)()()(1)()()aa(,2,3,,)(,2,3,,)kkkijijikkjkkkiiikklaijnbblbijn如此进行1n次消元，即1,2,,1kn，最后得同解方程组()()nnAxb其中，(1)(1)(1)(1)(1)11121311(2)(2)(2)(2)()()k()222322()()nnnknnnnnnaaaabaaabAxbAbab（）假若()0nnna，对式（3.1.3）进行回代计算得方程组的解()()()()1/()/(1,2,1)nnnnnniiiiiijjiijixbaxbaxain上述全过程称为Gauss顺序消元法。1.2Gauss主元素消元法Gauss顺序消元法每一步总有一个要求()0kkka，否则算法就失败。然而从方程组的理论，只要det()0A，则方程组解存在且唯一。由此可知，Gauss顺序消元法可执行条件苛刻。不仅如此，从数值计算角度，当()0kkka，但()kkka很小时。用()kkka做分母运算，会引起误差界的放大，数值不稳定现象产生，严重时导致解失真。为克服这一缺点，常采用选主元的消元法。Gauss列主元消元法主要依据线性方程组任意交换两个方程的次序，方程组的同解性不变，且解的分量次序也不变。于是，第k步在顺序消元法进行之前，从()kA的第k列元素()kkka，()1kkka，，()knka中选取绝对值最大者，并记录所在行，即()()maxkkkikikkinaa，记kli如果lk，则交换矩阵()()kkAb的第k行与第l行所有对应的元素，然后，再进行第k步顺序消元法。1.3Gauss全主元消元法Gauss全主元消元法是在Gauss列主元消元法的基础上进一步改进，目的是更好的提高数值稳定性和解的精度。对那些矩阵性态不太理想的方程组能给出最大可能性的满意解。全主元算法的第k步是从()(,1,,)kijaikkn中选取绝对值最大者作为主元，并记录主元所在列与所在行，即()(),maxkkkkijijkijnaa，记,kklitj如果lk，则交换矩阵()()kkAb的第k行与第l行所有对应的元素；如果tk，则交换矩阵()()kkAb的第k行与第t列所有对应的元素。值得注意的是，进行行元素的交换，解的分量次序不变，但是进行列元素交换时，解的分量次序将有所改变。1.4GaussSeidel迭代法将方程组Axb改写成等价方程组11211111111122212222222212000nnnnnnnnnnnnnaabaaaxxabaxxaaaxxbaaaaa(1.4.1)矩阵形式为JxGxf(1.4.2)其中1121111221222212000nnJnnnnnnaaaaaaaaGaaaa，111222nnnbabafba(1.4.3)任取初始解向量(0)(0)(0)(0)12(,,,),Tnxxxx带入式(1.4.2)右端得迭代格式(1)()(0,1,2,)kkJxGxfk其分量形式为1()1()/(1,2,,)(0,1,2,)nkkiiijjiijjixbaxaink(1.4.4)此式称为Jacobi迭代法，简称J法。将上式修改为11(1)()11()/(1,2,,)inkkkiiijjijjiijjixbaxaxain则为GaussSeidel迭代法。二、计算机的基本配置版本：7win旗舰版安装内存()RAM:2GB处理器：()IntelR()CoreTM3iCPUM3502.27GHz系统类型：32位操作系统三、参数设置及计算结果在本文中选取一个六元方程组12345612345612345612345612345612345623456442345787345469792123454653444122445894523433415235677435645xxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxx这是本人随机写的一个方程组，本文将以之为例，用摘要里说的方法来解方程组。3.1Gauss全主元消元法3.1.1算法流程图：1)输入(,1,2,,),(1,2,,)ijiaijnbin；2)()(1,2,,)ziiin；3)对1,2,,1kn做a)选主元：,maxkkijijkijnaa，并记,kklitj，即max;;;kkaalktk对1,,ikn做对1,,jkn做ifmaxijaathenmax,,ijaalitjb)ifmax0a，输出奇异标志，停机；c)iflkthen(,1,,);kjljklaajkknbb；iftkthen(1,2,,);()()ikitaainzkzt；d)对1,,ikn做i./ikikikkkalaa；ii.iiikkbblbiii.对1,,jkn做ijijikkjaala；4)if0nnathen输出奇异标志，并停机else做a)/nnnnbba；b)1()/(1,,1);niiijjiijibbabainc)(())(1,2,,)ixzibin；5)输出解(1,2,,)ixin3.1.2计算结果及分析（程序见附录）：123x=18.5115x7.2791x=-12.0967456x=-4.6715x-5.9736x=15.9624系数矩阵为：123456234578454697a=123454653444244589452343152356774356447392b=123445构成增广矩阵为123456442345787345469792123454653444122445894523433415235677435645A经过Gauss全主元消元法消元后矩阵a变为89.0045.000045.000023.000043.000024.000048.6854-5.314628.528128.9438-0.10110010.8117-2.0441-4.5008-2.48350006.28a064.12303.713300001.9065-1.414800000-0.2568消元后的b变为34.00023.6067-26.907797.03424.2430-4.7529b最终将解得的x带入原方程组，计算出误差矩阵errorAXb下面给出误差矩阵120.01420.02840.0284100.11370.11370.1137error由误差矩阵可知误差已接近于0，可知：这种算法的具有数值稳定，精度高的优点，是比较好的算法，但是缺点是在选主元的过程中，消耗了较多的时间。该实例的运行时间为0.048152seconds.3.2GaussSeidel迭代法3.2.1GaussSeidel迭代法算法流程图：1)输入：(,1,2,,),(1,2,,),;ijiaijnbin2)输入：初始解向量(1,2,,)ixin；3)对1,2,,in做a)1()/niiijjiijjixbaxa；b)iiieyx；c)iixy；4)max{}iiinfethen输出：ix(1,2,,)in，停机else返回第三步。3.2.2参数的设置：本题初始解：(0)(0,0,0,0,0,0)Tx;误差精度0.0001;3.2.3计算结果及分析（程序见附录）：针对这个线性方程组，用迭代法无法求解，因为这个系数矩阵的谱半径不小于1，在本题中谱半径为161.1222远大于1。这个问题的出现从侧面反映了迭代法的适应面会很小，因为现实生活中的很多矩阵都不是会符合谱半径小于1这个条件的，虽然这种方法迭代速度回很快，但是这个谱半径要求小于1弊端是其致命的弱点。四、计算结果分析及遇到的问题4.1计算结果的分析根据误差矩阵来看，此方法解决此类问题具有很好的效果，误差数量级在1210上面，误差是非常小的，几乎可以