1、3.2 Jacobi3.2 Jacobi迭代法和迭代法和Gauss-SeidelGauss-Seidel迭代法迭代法本节主要内容本节主要内容 Jacobi Jacobi迭代法迭代法 Gauss-Seidel Gauss-Seidel迭代法迭代法先看下面一个简单的例子.考虑线性方程组例例3.2.13.2.13.2.1 Jacobi3.2.1 Jacobi迭代法迭代法1231010911027,02108xxx 其精确解为T*1,1,1.x 【解解】把线性方程组改写为121232310 91027 2108xxxxxxx1221332191010117105101455xxxxxxx取=0,0,0
2、，由上式可构造迭代公式如下0Tx进一步改写为1121213132191010117,0,1,2,105101455kkkkkkkxxxxxkxx1310.kkxx直到计算结果见下表：1kxT*1,1,1.x 从计算结果看出，向量序列收敛到方程组的精确解kk000030.99500.98500.990010.90000.70000.800040.99850.99750.997020.97000.95000.940050.99980.99920.9995k1xk2xk3xk1xk2xk3x2kx3kx3kx1kx2kxJacobiJacobi迭迭代代法法.以上迭代法称为Carl Gustav Ja
3、cobi,1804-18511804.JacobiAbelJacobi卡尔雅可比（）是一位德国数学家.年生于波茨坦是历史上最伟大的数学家之一，他在数学方面的突出成就是和挪威数学家 相互独立地奠定了椭圆函数论的基础.的工作还包括代数学、变分法、复变函数论、微分方程和数学史等方面.Jacobi 在数值计算方面的主要贡献是提出求解线性方程组的迭代法以及求解矩阵的特征值和特征向量的方法等.卡尔卡尔.雅可比雅可比下面考虑一般的情形：(3.2.1)方程组的分量形式为0,iia 因所以有,Axb其中A是n阶非奇异矩阵.且其主对角元素0 (1,2,)iiain 0 iia (1,2,).in1,1,2,.ni
4、jijja xbin1,1,1,2,.niijjjj iiiixa xbina(3.2.1)(3.2.2)其中从而得到Jacobi迭代法的分量形式：11,1,1,2,.nkkiijjjj iiiixa xbina下面推导Jacobi迭代法的矩阵形式：把系数矩阵A分解成三部分：,ADLU1122diag(,),nnDaaa(3.2.4)(3.2.3)(3.2.5)(3.2.1)于是方程组改写为即211,11,212,100,00nnnnn naLaaaaa121,112,121,0000nnnnnnaaaaaUa()DxLU xb11().xDLU xD b11()JDLUID A任取向量 ，则
5、Jacobi迭代法可写成如下的矩阵形式：111(),0,1,2,kkxDLU xD bk111(),BDLUID AJfD b若记，则有1kkxJ xf称矩阵J 为Jacobi迭代法的迭代矩阵迭代矩阵.(3.2.8)0 x3.2.3.2.2 2 算法与程序算法与程序2.12.2.当时，执 行 步 骤kN步步 骤骤 2 2,1Nk精 度最 大 迭 代 次 数.令.算法算法3.1 Jacobi迭代法说明：为简单起见，假定系数矩阵A非奇异，且 ，且假设Jacobi迭代法收敛.0 (1,2,)iiain步骤步骤1 输入系数矩阵A，右端向量b，以及初始向量 0,x1,2,对，计 算in步步 骤骤 2 2
6、.1 111,1.nkkiijjijj iiixa xba111:,:1 若，则 算 法 停 止，输 出方 程 组 的 近 似 解 x；否 则，令.kkkkkxxxxkk 步步 骤骤 2 2.2 2步步 骤骤 3 3输 出 信 息“算 法 超 出 最 大 迭 代 次 数！”，算 法终 止.算法算法3.3.1 1的的 M Matlab atlab 程序程序M atlab 程程 序序 如如 下下：%Jacobi.mfunction x=Jacobi(A,b,x0,eps,N)%功能：用Jacobi迭代法解n 阶线性方程组 Ax=b n=length(b);x=ones(n,1);k=0;while
7、 k=N for i=1:n%步骤 2 x(i)=(b(i)-A(i,1:i-1,i+1:n)*x0(1:i-1,i+1:n)/A(i,i);end k=k+1 If norm(x-x0,inf)N warning(算法超出最大迭代次数！)；else disp(迭代次数=,num2str(k)x end例例3.2.23.2.2用Jacobi迭代法程序Jacobi.m求解线性方程组：12341012061111325.21 10111031815xxxx编写M文件调用函数Jacobi.m，并运行【解解】clc clear allformat longA=10,-1,2,0;-1,11,-1,3;
8、2,-1,10,-1;0,3,-1,8;b=6;25;-11;15;x0=0;0;0;0;eps=1e-3;N=300;x=Jacobi(A,b,x0,eps,N);计算结果为 迭代次数=10 x=1.000118598691415 1.999767947010035 -0.999828142874476 0.9997859784600503.2.3.2.3 3 Gauss-Seidel Gauss-Seidel迭代法迭代法Jacobi(3.2.2)将迭代法的迭代公式改写为11111,1,2,.inkkkiijjijjijj iiixa xa xbina (3.2.9)11(1)kkixixi
9、i 注意到计算的第 个分量时前面的 个分量已经计算出来，如果计算第 个分量时，前个分量使用最新的值，则得到111111,1,2,.inkkkiijjijjijj iiixa xa xbina (3.2.10)GaussSeidel.上式即为迭代法的分量形式(1)i GaussSeidel 下面推导迭代法的矩阵形式：(3.2.3),(3.2.1)据将方程组改写为(),DLU xb ,DxLxUxb 即(3.2.10)根据式，取迭代公式为11,kkkDxLxUxb111 ()(),kkxDLUxDLb 即GaussSeidel.这就是迭代法的矩阵形式Gauss-SeidelGauss-Seidel
10、迭代法的矩阵形式迭代法的矩阵形式 111()(),()若记，则有GDLUIDLAfDLb 1,kkxGxf GaussSeidel.G 迭迭代代矩矩阵阵称 为迭代法的Gauss-SeidelGauss-Seidel迭代法的迭代矩阵迭代法的迭代矩阵 11()()GDLUIDLA 例例 3.2.3 3.2.3 11211213113219,1010117,1051014,55kkkkkkkxxxxxxxT13010.0,0,03.2.kkxxx直至取初始向量，迭代结果见表GaussSeidel3.2.1 用迭代法解例的线性方程组，其迭代公式为表表3.23.2k1xk2xk3xk1xk2xk3xkk
11、000030.995700.997850.9995710.900000.700000.8000040.999790.999890.9999820.970000.957000.9914051.000000.999991.00000 T*53*3.2101,1,1.GaussSeidelJacobi.从表可知，其中是方程组的精确解 对于本例子，迭代法比迭代法收敛得快xxx 3.2.4 Gauss-Seidel3.2.4 Gauss-Seidel迭代法的步骤与程序迭代法的步骤与程序GaussSeidel下面给出迭代法的具体步骤：G Ga au us ss s-S Se ei id d算算法法e e3
12、 3.2 2 l l迭迭代代法法0 (1,2,)GaussSeidel.iiAain说明：为简单起见，假定系数矩阵 非奇异，且，且假设迭代法收敛0 ,.0.步步骤骤1 1 输入系数矩阵，右端向量，以及初始向量，精度最大迭代次数令AbxNk ,2.步步骤骤 2 2 当 时 执行步骤2.1-2.kN1,2,步步骤骤 2.12.1对，计算in111111:.inkkkiijjijjijj iiixa xa xba !.步步骤骤 3 3 输出信息 超出最大迭代次数，算法终止111:,:1 步步 骤骤 2 2.2 2 若，则算法停止，输出方程组的近似解 x；否则，令.kkkkkxxxxkk算法算法3.2
13、3.2的的M Matlabatlab程序程序%Gauss_Seidel.mfunction x=Gauss_Seidel A,b,x0,eps,N%GaussSeidel.nAxb功能：用迭代法解 阶线性方程组 nlength b;xones n,1;k0;while kNfor i1:n x iA i,1:i1*x 1:i1A i,i1:n*x0i1:nb i /A i,i;.endkk1;if norm xx0,infepsbreakend，；x0 x;endif kNWarning()算法超出最大迭代次数！；elsedisp(=,num2str(k)迭代次数endx 例例3.2.4 3.2.4 GaussSeidelGauss_Seidel.m3.2.2.用迭代法程序求解例的线性方程组 MGauss_Seidel.m,.编写如下文件调用函数并运行解解clcclear allformat longA10,1,2,0;1,11,1,3;2,1,10,1;0,3,1,8;b6;25;11;15;x00;0;0;0;eps1e3;N300;xGauss_Seidel A,b,x0,eps,N;
