解方程组基本思想.docx

1、解方程组基本思想四：基本方法基本思路将在解题的过程中得到体现。1（求线性方程组的唯一解或特解），这类问题的求法分为两类：一类主要用于解低阶稠密矩阵 直接法；一类是解大型稀疏矩阵 迭代法。1.1利用矩阵除法求线性方程组的特解（或一个解）方程：AX=b，解法：X=Ab，(注意此处不是/)例1-1 求方程组 的解。解: A = ; = ；b=(1,0,0,0,1)由于rank(A)=5,rank( )=5 %求秩，此为R（A）=R（ ）=n的情形，有唯一解。X= Ab %求解 X =(2.2662, -1.7218, 1.0571,-0.5940, 0.3188) 或用函数rref求解，sv=rre

2、f(A:b);所得sv的最后一列即为所要求的解。12 利用矩阵的LU、QR和cholesky分解求方程组的解 这三种分解，在求解大型方程组时很有用。其优点是运算速度快、可以节省磁盘空间、节省内存。I) LU分解又称Gauss消去分解，可把任意方阵分解为下三角矩阵的基本变换形式（行交换）和上三角矩阵的乘积。即A=LU，L为下三角阵，U为上三角阵。则：A*X=b 变成L*U*X=b所以X=U(Lb) 这样可以大大提高运算速度。命令 L，U=lu (A)在matlab中可以编如下通用m 文件：在Matlab中建立M文件如下% exp1.mA;b;L，U=lu (A);X=U(Lb)II）Choles

3、ky分解若A为对称正定矩阵，则Cholesky分解可将矩阵A分解成上三角矩阵和其转置的乘积，即： 其中R为上三角阵。方程 A*X=b 变成 所以在Matlab中建立M文件如下% exp2.mA;b;R，R=chol(A);X=R(Rb)III）QR分解对于任何长方矩阵A，都可以进行QR分解，其中Q为正交矩阵，R为上三角矩阵的初等变换形式，即：A=QR方程 A*X=b 变形成 QRX=b所以 X=R(Qb)上例中 Q, R=qr(A)X=R(QB)在Matlab中建立M文件如下% exp3.mA;b;Q，R=qr(A);X=R(Qb)2求线性齐次方程组的通解(A*X=0)在Matlab中，函数n

4、ull用来求解零空间，即满足A•X=0的解空间，实际上是求出解空间的一组基（基础解系）。在Matlab中建立M文件如下% exp4.mformat rat %指定有理式格式输出A;b=0;r=rank(A);bs=null(A，r); %一组基含(n-r)个列向量% k ,k ,k % X= k *bs(:,1)+ k *bs(:,2)+ k *bs(:,n-r) 方程组的通解pretty(X) %让通解表达式更加精美3 求非齐次线性方程组的通解（A*X=b）非齐次线性方程组需要先判断方程组是否有解，若有解，再去求通解。因此，步骤为：第一步：判断AX=b是否有解，(利用基本思路的第

5、一条)若有解则进行第二步第二步：求AX=b的一个特解第三步：求AX=0的通解第四步：AX=b的通解为： AX=0的通解加上AX=b的一个特解。在Matlab中建立M文件如下% exp4.mclear allA;b; %输入矩阵A,bm,n=size(A);R=rank(A);B=A b;Rr=rank(B);format rat if R=Rr&R=n % n为未知数的个数，判断是否有唯一解x=Ab;elseif R=Rr&R |aij| i=1,2,n, j=1,j i则称方阵A是严格(行)对角占优的.7.收敛定理 对任意初始向量x(0)及任意右端向量 g，由迭代x(k+1) =B x(k)

6、 +g产生的迭代向量序列x(k)收敛的充要条件是谱半径 （B）18.收敛判别条件判别条件1： 若|B|1, 则迭代x(k+1) =B x(k) +g 对任何初始向量x(0)都收敛. 判别条件2： 如果A为严格对角占优阵，则其 Jacobi迭代和Seidel迭代对任何初始向量x(0)都收敛。 判别条件3： 如果A为对称正定阵，则其 Seidel迭代对任何初始向量x(0)都收敛。9.迭代法的误差估计 若|B|1，则对迭代格式 x(k+1) =B x(k) +g 有3.3 程序中Mathematica语句解释a*matrix 数a与矩阵matrix相乘matrix1+matrix2 矩阵matrix

7、1和矩阵matrix2相加（注意矩阵的大小相同）matrix1.matrix2 矩阵matrix1和矩阵matrix2相乘(注意矩阵乘法的规则) Transposematrix 求矩阵matrix转置Inversematrix 求矩阵(方阵) matrix 的逆 DiagonalMatrixlist使用列表list中的元素生成一个对角矩阵.IdentityMatrixn 生成n阶单位矩阵Maxx求向量x中元素的最大值3.4方法、程序、实验 解线性方程组的迭代法是将线性方程组 Ax=b 化为等价线性方程组 x=Bx+f 再由矩阵迭代格式 x(k+1)=Bx (k)+f构造向量序列x(k)来求线性

8、方程组解的。如果得出的向量序列x(k)收敛至某个向量x*，则可得该向量x*就是所求方程组 Ax=b 的准确解.线性方程组的迭代法主要有Jocobi迭代法、Seidel迭代法和超松弛(Sor)迭代法。1. Jocobi迭代法1) Jocobi迭代法的构造过程 假设aii 0，依次在第i个方程解出x i , i=1,2, ,n并令 cij = -aij /aii (i j) ，gi= bi /aii 就得到如下Jocobi迭代格式：x1(k+1)= c12x2(k)+c13x3(k)+ +c1nxn(k)+g1 x2(k+1)=c21x1(k) +c23x3(k)+ +c2nxn(k)+g2 。

9、xn(k+1)=1x1(k) +2x2(k)+ +(n-1)xn-1(k) + gn若令则有Jocobi迭代的矩阵格式：x(k+1) = BJx(k) +gJBJ 称为Jocobi迭代矩阵。Jocobi迭代可以写成如下紧凑格式：在给定初始迭代向量x(0)后就可以进行Jocobi迭代求解了。2)Jacobi迭代算法1.输入变量个数n、初值向量x（0）、迭代精度eps、系数矩阵A、常数项b 和迭代最大次数nmax2 For i=1,2,n2.1 如果|aii|eps1,则输出“迭代失败”提示并终止3. Bj E-D-1A4. gj D-1b5.For k=1,2,nmax5.1 x Bj.x0+

10、gj5.2 如果|x-x0|eps ,输出迭代失败，终止。3)Jacobi 迭代法程序Cleara,b,x;nmax=500;n=Input“线性方程组阶数n=”;a=Input系数矩阵A=；b=Input常数项b=；x0=Input输入迭代初值向量x0；eps1=0.000001;eps=Input输入精度控制eps=；DoIfAbsai,ieps1,t1=1;Return,t1=0,i,1,n;Ift1=1, PrintJacobi迭代法失效, d=DiagonalMatrixTableai,i,i,1,n; d1=Inversed; bj=IdentityMatrixn-d1.a; gj

11、=d1.b; Do x=bj.x0+gj; err=MaxAbsx-x0; Printx=,x/N, i=,i, err=,err/N; IfNerr=eps,Print迭代失败 说明 本程序用于求线性方程组Ax=b的解。程序执行后，先通过键盘输入线性方程组阶数n、系数矩阵A、常数项b、迭代初值向量x0和输入精度控制eps，程序即可给出每次迭代的次数和对应的迭代向量序列x(k),其中最后输出的结果即为所求的根。如果迭代超出500次还没有求出满足精度的根则输出迭代失败提示，如果出现主对角线元素aii=0给出Jacobi迭代法失效提示。程序中变量说明x0：存放初始向量和迭代过程中的向量x（k）x:

12、 存放迭代过程中的向量x（k+1）nmax:存放迭代允许的最大次数err:存放误差|x-x0| t1:临时变量注：迭代最大次数可以修改为其他数字。4)例题与实验例1.用Jacobi 迭代法解如下线性方程组 5x1+2x2+x3= -12 -x1+4x2+2x3= 20 2x1-3x2+10x3= 3要求误差|x(k+1)-x(k)| 10-4，并用取不同初值的方法实验观察迭代收敛的情况。解：执行Jacobi迭代法程序后在输入的四个窗口中按提示分别输入:3、5, 2, 1, -1, 4, 2, 2, -3, 10、-12, 20, 3、0, 0, 0、0.0001每次输入后用鼠标点击窗口的“OK

13、”按扭，得如下输出结果： x=-2.4, 5., 0.3 i=1 err=5.x=-4.46, 4.25, 2.28 i=2 err=2.06x=-4.556, 2.745, 2.467 i=3 err=1.505x=-3.9914, 2.6275, 2.0347 i=4 err=0.5646x=-3.85794, 2.9848, 1.88653 i=5 err=0.3573x=-3.97123, 3.09225, 1.96703 i=6 err=0.113286x=-4.03031, 3.02368, 2.02192 i=7 err=0.0685705x=-4.01386, 2.98146,

14、 2.01316 i=8 err=0.042216x=-3.99522, 2.98995, 1.99721 i=9 err=0.0186374x=-3.99542, 3.00259, 1.99603 i=10 err=0.0126367x=-4.00024, 3.00313, 1.99986 i=11 err=0.0048186x=-4.00122, 3.00001, 2.00099 i=12 err=0.00312067x=-4.0002, 2.9992, 2.00025 i=13 err=0.00102319x=-3.99973, 2.99983, 1.9998 i=14 err=0.00

15、0625697x=-3.99989, 3.00017, 1.99989 i=15 err=0.00034136x=-4.00005, 3.00008, 2.00003 i=16 err=0.000155236x=-4.00004, 2.99997, 2.00003 i=17 err=0.000106099x=-4., 2.99997, 2. i=18 err=0.0000414468此结果说明迭代18次，求得误差为err=0.0000414468的近似解，最后显示的近似解向量为x=-4., 2.99997, 2.，它表示所求解为 x1=-4，x2=2.99997，x3=2 。本题的准确解为 x

16、1=-4，x2=3，x3=2 。如果将如上输入的初值改为21, -18, 30，执行Jacobi迭代法程序后得如下输出结果： x=-1.2, -4.75, -9.3 i=1 err=39.3x=1.36, 9.35, -0.885 i=2 err=14.1x=-5.963, 5.7825, 2.833 i=3 err=7.323x=-5.2796, 2.09275, 3.22735 i=4 err=3.68975x=-3.88257, 2.06642, 1.98375 i=5 err=1.39703x=-3.62332, 3.03749, 1.69644 i=6 err=0.97106x=-3

17、.95428, 3.24595, 1.93591 i=7 err=0.330963x=-4.08556, 3.04347, 2.06464 i=8 err=0.202475x=-4.03032, 2.94629, 2.03015 i=9 err=0.0971858x=-3.98455, 2.97734, 1.98995 i=10 err=0.0457716x=-3.98893, 3.00889, 1.99011 i=11 err=0.0315451x=-4.00158, 3.00771, 2.00045 i=12 err=0.0126504x=-4.00318, 2.99938, 2.0026

18、3 i=13 err=0.00833246x=-4.00028, 2.99789, 2.00045 i=14 err=0.00289753x=-3.99925, 2.99971, 1.99942 i=15 err=0.0018145x=-3.99977, 3.00048, 1.99976 i=16 err=0.000770763x=-4.00014, 3.00018, 2.0001 i=17 err=0.000375926x=-4.00009, 2.99992, 2.00008 i=18 err=0.000261658x=-3.99998, 2.99994, 1.99999 i=19 err=

19、0.000107579x=-3.99997, 3.00001, 1.99998 i=20 err=0.0000714045从计算结果可以看到虽然所取的初值不同，但迭代总是收敛相同的结果。考察本题系数矩阵可以发现它是严格对角占优矩阵，由收敛判别法，Jacobi迭代对任意初值都是收敛的，我们通过实际计算验证了这个结果。例2.通过Jacobi 迭代序列观察用Jacobi 迭代法解如下线性方程组 x1+2x2+2x3= -12 -x1+4x2+x3= 20 2x1-3x2+x3= 3的收敛性。解：执行Jacobi迭代法程序后在输入的四个窗口中按提示分别输入:3、1, 2, 2, -1, 4, 1, 2

20、, -3, 1、-12, 20, 3、0, 0, 0、0.0001每次输入后用鼠标点击窗口的“OK”按扭，得如下输出结果： x=-12., 5., 3. i=1 err=12.x=-28., 1.25, 42. i=2 err=39.x=-98.5, -12.5, 62.75 i=3 err=70.5x=-112.5, -35.3125, 162.5 i=4 err=99.75.246.719, -120.176, 696.719 i=8 err=763.5x=-1165.09, -107.5, -850.965 i=9 err=1547.68x=1904.93, -73.5303, 2010

21、.67 i=10 err=3070.02x=-3886.28, -21.4355, -4027.45 i=11 err=6038.12x=8085.77, 40.2917, 7711.26 i=12 err=11972.1x=-15515.1, 98.6279, -16047.7 i=13 err=23758.9x=31886.1, 138.141, 31329.1 i=14 err=47401.2x=-62946.5, 144.247, -63354.7 i=15 err=94832.5x=126409., 107.068, 126329. i=16 err=189683.x=-252883

22、., 25.0775, -252494. i=17 err=379292.x=504925., -92.4305, 505845. i=18 err=758339.通过观察迭代过程中的误差是不断变大的特点可以知道本题的Jacobi 迭代序列是不收敛的，因此，本题线性方程组不能用Jacobi 迭代法求解。这个实验说明并不是每个线性方程组都能用Jacobi 迭代法求解。2. Seidel迭代1) Seidel迭代的构造过程 为了加快收敛速度，同时节省计算机的内存，对Jocobi迭代作如下的改进：每算出一个分量的近似值，立即用到下一个分量的计算中去，即用迭代格式： x1(k+1)= c12x2(k)

23、+c13x3(k)+ +c1nxn(k)+g1 x2(k+1)=c21x1(k+1) +c23x3(k)+ +c2nxn(k)+g2 。 xn(k+1)=1x1(k+1) +2x2(k+1)+ +(n-1)xn-1(k+1) + gn如上迭代可以写成如下紧凑格式：这样所得的迭代法就称为Seidel迭代法。 利用Ax=b 及A=L+D+U,其中D为对角矩阵,L,U分别为严格下，上三角矩阵则有Seidel迭代法的矩阵形式为 x(k+1)= Bs x(k)+gs 其中：Bs =-(L+D)-1U，g s=D-1bSeidel迭代矩阵为 Bs =-(L+D)-1U。 在给定初始迭代向量x(0)后就可以

24、进行Seidel迭代求解了。 Jacobi迭代和Seidel迭代格式可表述为统一形式:2) Seidel迭代算法1.输入变量个数n、初值向量x（0）、迭代精度eps、系数矩阵A、常数项b 和迭代最大次数nmax2. For i=1,2,n2.1 如果|aii|eps1,则输出“迭代失败”提示并终止3.For i=1,2,nmax 3.1 For i=1,2,n3.2 如果|x(k+1)-x(k)| eps ,输出解向量x，终止；否则x（k） x(k+1)4. 如果|x-x0| eps ,输出迭代失败，终止。3) Seidel 迭代法程序Cleara,b,x;nmax=500;n=Input“线

25、性方程组阶数n=”;a=Input系数矩阵A=；b=Input常数项b=；x0=Input输入迭代初值向量x0；eps1=0.000001;eps=Input输入精度控制eps=； x=x0; DoIfAbsai,ieps1,t1=1;Return,t1=0,i,1,n;Ift1=1, PrintSeidel迭代法失效, Do Dou1=Sumai,j*xj,j,1,i-1; u2=Sumai,j*x0j,j,i+1,n; xi=(bi-u1-u2)/ai,i, i,1,n; err=MaxAbsx-x0; Printx=,x/N, k=,k, err=,err/N; IfNerr=eps,Print迭代失败 说明 本程序用于求线性方程组Ax=b的解。程序执行后，先通过键盘输入线性方程组阶数n、系数矩阵A、常数项b、迭代初值向量x0和输入精度控制eps，程序即可给出每次迭代的次数和对应的迭代向量序列x(K),其中最后输出的结果即为所求的根。如果迭代超出500次还没有求出满足精度的根则输出迭代失败提示，如果出现主对角线元素aii=0给出Jacobi迭代法失效提示。程序中变量说明x0：存放初始向量和迭