Jacobi迭代法和Gauss-Seidel迭代法Matlab程序
- 格式:docx
- 大小:13.83 KB
- 文档页数:6
Jacobi与Gauss-Seidel迭代的比较及r算法的MATLAB实现胡志成【摘要】探讨比较了Jacobi迭代和Gauss-Seidel迭代,尤其是它们的可并行性.给出了它们在MATLAB中的一个实现.通过算例说明,算法的并行化实现可以大幅度地提升Jacobi迭代的效率.【期刊名称】《高师理科学刊》【年(卷),期】2018(038)003【总页数】4页(P59-61,84)【关键词】Jacobi迭代;Gauss-Seidel迭代;并行性【作者】胡志成【作者单位】南京航空航天大学理学院,江苏南京 210016【正文语种】中文【中图分类】O241.6;G642.0Jacobi迭代和Gauss-Seidel迭代是解线性方程组的2个经典的迭代法,简单又易于实现,是大学计算方法相关课程必讲的知识点.然而教材[1-3]往往只注重理论而轻视算法的程序实现.就这2个迭代法而言,表现为着重讨论迭代法的收敛性和收敛速度,而甚少给出或者仅给出一个粗糙低效的算法实现.至于Jacobi迭代的并行性,探讨的更是少之又少.由于课时的限制,在工科计算方法教学中更是如此.这一方面是由于一般程序语言关于算法的并行化实现并不那么容易,另一方面则是因为在课堂上短时间内举例说明算法并行化实现带来的效率提升较为困难.然而在课堂上适当地给出算法的实现和具体算例,显然有助于学生对算法本身的理解.本文首先回顾了Jacobi迭代和Gauss-Seidel迭代的基本格式和传统比较,然后探讨了它们的可并行性,并给出了在MATLAB中的一个实现.文中给出的算例不仅可以在课堂的短时间内展示结果,也能说明算法的并行化实现可以大幅度地提升Jacobi迭代的效率.设是阶非奇异方阵,为维列向量.若的主对角元,,则线性方程组可以用Jacobi 迭代和Gauss-Seidel迭代求解.其迭代格式分别为若再将表示成,其中:;;,则Jacobi迭代可用矩阵表示为相应地,Gauss-Seidel迭代的矩阵形式为关于Jacobi迭代和Gauss-Seidel迭代的分析和比较,在众多的文献中都有提及[1,4-6],主要归结为:(1)收敛性.当且仅当各自的迭代矩阵(分别是和)的谱半径小于1时,该迭代法收敛.若是严格对角占优阵,则2个迭代法都收敛.但一般地, Jacobi迭代收敛不能保证同一问题的Gauss-Seidel迭代也收敛,反之亦然.(2)收敛速度.当Gauss-Seidel迭代和Jacobi迭代两者都收敛时,多数情况下Gauss-Seidel迭代比Jacobi迭代收敛得更快.对于一些模型问题,Gauss-Seidel迭代的收敛速度甚至是Jacobi迭代的2倍以上[5]58.也就是说,从同一初值出发,为得到满足精度要求的解,Gauss-Seidel迭代所需的步数可能比Jacobi迭代所需步数的一半还要少.(3)每步计算量.由于矩阵的求逆并不那么容易,具体的算法实现多采用分量形式,即式(1)和式(2).比较式(1)和式(2)可知,2个迭代法每步的计算量相等.(4)存储量.Gauss-Seidel迭代只需1个向量来存储或的分量,而Jacobi迭代需要2个向量来分别存储和.比较可知,在传统的算法实现中,Gauss-Seidel迭代比Jacobi迭代往往更具优势.随着多核处理器的迅猛发展,并行计算渐次成为提升计算效率的重要技术手段,因而一个算法的可并行性也成为了衡量其优劣的重要指标[7-9].进一步比较迭代法(1)和(2),易知Jacobi迭代中,各分量的计算是相互独立的,因而具有高度的并行性;而Gauss-Seidel迭代中,各分量是依次顺序计算的,因而本质上是一个串行算法.这就意味着,利用算法的可并行性,如果实现得当,即使收敛速度较慢,Jacobi迭代仍然有可能在更短的时间内得到满足精度要求的数值解.为了加深学生对此的理解,并让学生体会到算法实现,特别是并行化实现的重要性,有必要给出适当的实例和程序实现.MATLAB是一种面向科学与工程计算的数学软件,因其强大的数值计算和图形绘制能力、广泛的工具包支持以及使用的便利性,深受教学和科研工作者的喜爱,是工科大学生需要掌握的重要软件之一[10].借助其灵活的矩阵向量运算,给出Jacobi迭代(1)在MATLAB环境中的一个实现(见算法1)和Gauss-Seidel迭代(2)的一个MATLAB实现(见算法2).可以看到,算法1完全避免了Jacobi迭代各分量之间的for循环,而算法2则无法避免Gauss-Seidel迭代各分量之间的for循环,这是导致2个算法执行效率巨大差异的直接原因.算法1(Jacobi 迭代):function [x,k,res] = Jacobi(A,b,x0,tol,kmax )k = 0; res = 1;while(k <= kmax)& res > tol)x =(b - A * x0 + diag(A).* x0 )./ diag(A);k = k+1;res = max(abs(b - A * x ));x0 = x;end算法2( Gauss-Seidel 迭代):function [x,k,res] = GaussSeidel(A,b,x0,tol,kmax )k = 0;res = 1;x = x0;while (k <= kmax)&(res > tol)for i = 1:length(x)x(i)=(b(i)- A(i,1:end)* x + A(i,i)* x(i))/ A(i,i);endk = k+1;res = max(abs(b - A * x));end以零向量为初值,以残量的无穷范数为收敛准则,考察比较这2个算法求解线性方程组的数值结果,包括迭代步数、每步耗时和总耗时.为了使测试算例具有一定的普遍意义,且能在课堂上较短的时间内展示比较结果.用算法3来随机生成阶线性方程组的系数矩阵.算法3(随机生成一个严格对角占优矩阵):function A = randSDDMatrix(n)nzmax = randi([n*5,round(n*10)],1);i = [randi([1,n],nzmax,1)];j = [randi([1,n],nzmax,1)];s = rand(nzmax,1);A = sparse(i,j,s,n,n);x = randi([1,1e1],n,1)+ sum( abs(A),2);ii = 1:n;A = sparse(ii,ii,x,n,n)+ A;该矩阵是稀疏的,并且是严格对角占优的,可以确保阶数较大时,2个迭代法仍能快速收敛.当然,右端项也是随机生成的.对每一个阶数,都测试100组和,并取其平均数据(见表1).由表1可以看出,虽然Jacobi迭代(1)和Gauss-Seidel迭代(2)的每步计算量相等,但是算法1的每步计算效率显著优于算法2.这最终导致:虽然Jacobi 迭代的总迭代步数是Gauss-Seidel迭代的2倍以上,但是总耗时却远少于Gauss-Seidel迭代.而且随着矩阵阶数增加,在本文的算法实现下,这个优势将进一步扩大.本文比较分析了Jacobi迭代和Gauss-Seidel迭代,强调说明Jacobi迭代是一个高度并行的算法,而Gauss-Seidel迭代则是一个串行算法.因而如果算法实现得当,即使收敛速度更慢的时候,Jacobi迭代仍然可能比Gauss-Seidel迭代更快地给出数值解.在教学实践中,课堂上展示这些MATLAB实现和数值算例,明显地加深了学生对这2个迭代法的理解,也使学生意识到了算法实现的重要性,获得了较好的教学效果.【相关文献】[1] 倪勤,王正盛,刘皞.数值计算方法[M].北京:高等教育出版社,2012[2] 李庆杨,王能超,易大义.数值分析[M].5版.北京:清华大学出版社,2008[3] 李有法,李晓勤.数值计算方法[M].2版.北京:高等教育出版社,2005[4] 杜衡吉,徐昆良.Jacobi和Gauss-Seidel迭代法求解线性方程组的分析及应用[J].曲靖师范学院学报,2011,30(3):46-50[5] 白红梅.Jacobi迭代法与Gauss-Seidel迭代法的收敛性比较分析[J].呼伦贝尔学院学报,2009,17(6):55-58[6] 包丽君.“数值分析”中实践教学的探讨[J].宁波广播电视大学学报,2008,6(2):91-93[7] 徐树方,高立,张平文.数值线性代数[M].北京:北京大学出版社,2000[8] 杨学军.并行计算六十年[J].计算机工程与科学,2012,34(8):1-10[9] 韩建勋.并行计算在矩阵运算中的应用[J].信息与电脑:理论版,2017(5):93-96[10] 韦慧,蒋利华.基于MATLAB的大学工科《计算方法》课程教学方法探讨[J].内江科技,2017(4):72-73。
Jacobi 迭代法与Gauss-Seidel迭代法算法比较目录1 引言 (1)1.1Jacobi迭代法 (2)1.2Gauss-Seidel迭代法 (2)1.3逐次超松弛(SOR)迭代法 (3)2算法分析 (3)3 结论 (5)4 附录程序 (5)参考文献 (8)Jacobi 迭代法与Gauss-Seidel 迭代法比较1 引言解线性方程组的方法分为直接法和迭代法,直接法是在没有舍入误差的假设下,能在预定的运算次数内求得精确解,而迭代法是构造一定的递推格式,产生逼近精确值的序列。
这两种方法各有优缺点,直接法普遍适用,但要求计算机有较大的存储量,迭代法要求的存储量较小,但必须在收敛性得以保证的情况下才能使用。
对于高阶方程组,如一些偏微分方程数值求解中出现的方程组,采用直接法计算代价比较高,迭代法则简单又实用,所以比较受工程人员青睐。
迭代法求解方程组就是构造一个无限的向量序列,使它的极限是方程组的解向量。
即使计算机过程是精确的,迭代法也不能通过有限次算术运算求得方程组的精确解,而只能逐步逼近它。
因此迭代法存在收敛性与精度控制的问题。
迭代法是常用于求解大型稀疏线性方程组(系数矩阵阶数较高且0元素较多),特别是某些偏微分方程离散化后得到的大型稀疏方程组的重要方法。
设n 元线性微分方程组b Ax = (1)的系数矩阵A 非奇异,右端向量0≠b ,因而方程组有唯一的非零解向量。
而对于这种线性方程组的近似解,前辈们发展研究了许多种有效的方法,有Jacobi 迭代法、Gauss —Seidel 迭代法,逐次超松弛迭代法(SOR 法),这几种迭代方法均属一阶线性定常迭代法,即若系数矩阵A 分解成两个矩阵N 和P 的差,即P N A -=;其中N 为可逆矩阵,线性方程组(1)化为:b x P N =-)(b Px Nx +=b N Px N x 11--+=可得到迭代方法的一般公式:d Gx xk k +=+))1(( (2)其中:P N G 1-=,b N d 1-=,对任取一向量)0(x 作为方程组的初始近似解,按递推公式产生一个向量序列)1(x ,)2(x ,...,)k x(,...,当k 足够大时,此序列就可以作为线性方程组的近似解。
MATLAB块雅克比迭代法求线性方程组Ax=b的解块高斯-赛德尔迭代法求线性方程组Ax=b的解function [x,N]= BJ(A,b,x0,d,eps,M) %块雅克比迭代法求线性方程组Ax=b的解if nargin==4eps= 1.0e-6;M = 10000;elseif nargin<4errorreturnelseif nargin ==5M = 10000; %参数的默认值endNS = size(A);n = NS(1,1);if(sum(d) ~= n)disp('分块错误!');return;endbnum = length(d);bs = ones(bnum,1);for i=1:(bnum-1)bs(i+1,1)=sum(d(1:i))+1;%获得对角线上每个分块矩阵元素索引的起始值endDB = zeros(n,n);for i=1:bnumendb = bs(i,1)+d(i,1)-1;DB(bs(i,1):endb,bs(i,1):endb)=A(bs(i,1):endb,bs(i,1):endb );%求A的对角分块矩阵endfor i=1:bnumendb = bs(i,1)+d(i,1)-1;invDB(bs(i,1):endb,bs(i,1):endb)=inv(DB(bs(i,1):endb,bs(i ,1): endb));%求A的对角分块矩阵的逆矩阵endN = 0;tol = 1;while tol>=epsx = invDB*(DB-A)*x0+invDB*b; %由于LB+DB=DB-AN = N+1; %迭代步数tol = norm(x-x0); %前后两步迭代结果的误差x0 = x;if(N>=M)disp('Warning: 迭代次数太多,可能不收敛!');return;endendfunction [x,N]= BGS(A,b,x0,d,eps,M) %块高斯-赛德尔迭代法求线性方程组Ax=b的解if nargin==4eps= 1.0e-6;M = 10000;elseif nargin<4errorreturnelseif nargin ==5M = 10000;endNS = size(A);n = NS(1,1);bnum = length(d);bs = ones(bnum,1);for i=1:(bnum-1)bs(i+1,1)=sum(d(1:i))+1;%获得对角线上每个分块矩阵元素索引的起始值endDB = zeros(n,n);for i=1:bnumendb = bs(i,1)+d(i,1)-1;DB(bs(i,1):endb,bs(i,1):endb)=A(bs(i,1):endb,bs(i,1):endb ); %求A的对角分块矩阵endLB = -tril(A-DB); %求A的下三角分块阵UB = -triu(A-DB); %求A的上三角分块阵N = 0;tol = 1;while tol>=epsinvDL = inv(DB-LB);x = invDL*UB*x0+invDL*b; %块迭代公式N = N+1;tol = norm(x-x0);x0 = x;if(N>=M)disp('Warning: 迭代次数太多,可能不收敛!');return;endend类别:matlab 编程 | | 添加到搜藏 | 分享到i贴吧 | 浏览(168) | 评论 (0)上一篇:MATLAB 共轭梯度法求线性方程组A...。
【题目】:Gauss-Seidel迭代法及Matlab代码实例【内容】:1. Gauss-Seidel迭代法介绍Gauss-Seidel迭代法是一种用于解线性方程组的数值方法,基于逐次逼近的思想,通过不断迭代逼近线性方程组的解。
该方法通常用于求解大型稀疏线性方程组,其收敛速度相对较快。
2. 迭代公式推导假设有如下线性方程组:$$Ax=b$$其中A为系数矩阵,b为常数向量,x为未知向量。
Gauss-Seidel迭代法的迭代公式为:$$x^{(k+1)}=(D+L)^{-1}(b- Ux^{(k)})$$其中,D为A的对角矩阵,L为A的严格下三角矩阵,U为A的严格上三角矩阵,k为迭代次数。
3. Matlab代码实现下面给出Gauss-Seidel迭代法的Matlab代码实例:```matlabfunction [x, k] = gaussSeidel(A, b, x0, tol, maxIter)A: 系数矩阵b: 常数向量x0: 初始解向量tol: 容差maxIter: 最大迭代次数x: 解向量k: 迭代次数n = length(b);x = x0;k = 0;while k < maxIterx_old = x;for i = 1:nx(i) = (b(i) - A(i,1:i-1)*x(1:i-1) - A(i,i+1:n)*x_old(i+1:n)) / A(i,i); endif norm(x - x_old, inf) < tolreturnendk = k + 1;enddisp('迭代次数达到最大值,未达到容差要求'); end```4. 应用实例假设有如下线性方程组:$$\begin{cases}2x_1 - x_2 + x_3 = 5\\-x_1 + 2x_2 - x_3 = -2\\x_1 - x_2 + 2x_3 = 6\end{cases}$$系数矩阵A为:$$\begin{bmatrix}2 -1 1\\-1 2 -1\\1 -1 2\end{bmatrix}$$常数向量b为:$$\begin{bmatrix}5\\-2\\6\end{bmatrix}$$取初始解向量x0为:$$\begin{bmatrix}0\\0\\\end{bmatrix}$$容差tol为1e-6,最大迭代次数maxIter为100。
分别用 jacobi 迭代法和 gauss-seidel 迭代法,求解方程组【jacobi 迭代法和 gauss-seidel 迭代法分别应用于方程组的求解】1. 引言在数学领域中,方程组的求解一直是一个重要的课题。
为了解决复杂的线性方程组,人们提出了各种迭代方法,其中 jacobi 迭代法和gauss-seidel 迭代法是两种常见的方法。
本文将探讨这两种迭代方法在求解方程组中的应用。
2. jacobi 迭代法的原理和应用jacobi 迭代法是一种基于逐次逼近的迭代方法。
对于线性方程组AX=B,其中 A 是系数矩阵,X 是未知数向量,B 是已知向量。
我们可以通过以下公式进行逐次逼近:X(k+1) = D^(-1)*(B - (L+U)X(k))其中,D、L、U 分别是 A 的对角线、下三角和上三角矩阵。
jacobi 迭代法的优点在于易于理解和实现,但在收敛速度上较慢,需要进行多次迭代才能得到精确解。
在实际应用中,需要根据实际情况选择合适的迭代次数。
3. gauss-seidel 迭代法的原理和应用与 jacobi 迭代法类似,gauss-seidel 迭代法也是一种基于逐次逼近的迭代方法。
不同之处在于,gauss-seidel 迭代法在计算 X(k+1) 时利用了已经得到的 X(k) 的信息,即:X(k+1)_i = (B_i - Σ(A_ij*X(k+1)_j,j≠i))/A_ii这种方式使得 gauss-seidel 迭代法的收敛速度较快,通常比 jacobi 迭代法更快,尤其是对于对角占优的方程组。
4. 分别用 jacobi 迭代法和 gauss-seidel 迭代法求解方程组为了更具体地说明 jacobi 迭代法和 gauss-seidel 迭代法的应用,我们分别用这两种方法来求解以下方程组:2x1 + x2 = 9x1 + 3x2 = 11我们将该方程组写成矩阵形式 AX=B:|2 1| |x1| |9||1 3| * |x2| = |11|我们根据 jacobi 迭代法和 gauss-seidel 迭代法的原理,依次进行迭代计算,直到满足收敛条件。
Gauss-Seidel迭代法是一种用于解线性方程组的数值方法,特别适用于稀疏矩阵。
以下是一个使用Matlab实现Gauss-Seidel迭代法的简单示例代码:```matlabfunction [x, iteration] = gaussSeidel(A, b, tol, maxIter)% 输入参数:% A:系数矩阵% b:右侧常数向量% tol:迭代收敛容差% maxIter:最大迭代次数n = length(b);x = zeros(n, 1); % 初始化解向量iteration = 0;while iteration < maxIterx_new = x; % 存储前一次迭代的解for i = 1:n% 计算新的解x(i) = (b(i) - A(i,1:i-1)*x(1:i-1) - A(i,i+1:n)*x_new(i+1:n)) / A(i,i);end% 计算迭代误差error = norm(x - x_new, inf);if error < tolreturn;enditeration = iteration + 1;endwarning('Gauss-Seidel迭代未收敛到指定容差。
');end```使用这个函数时,您需要提供系数矩阵A、右侧常数向量b、迭代收敛容差tol和最大迭代次数maxIter作为输入参数。
函数将返回解向量x和迭代次数iteration。
示例用法:```matlabA = [4, -1, 0; -1, 4, -1; 0, -1, 3];b = [12; -1; 0];tol = 1e-6;maxIter = 1000;[x, iteration] = gaussSeidel(A, b, tol, maxIter);fprintf('解向量x = \n');disp(x);fprintf('迭代次数= %d\n', iteration);```这将求解线性方程组Ax = b,并返回解向量x以及迭代次数。
实验1:线性方程组的迭代解法1、实验环境MATLAB2009A2、实验目的和要求目的:利用Gauss-Seidel编程法求解方程组要求:代码能列出每一次迭代的中间值3、解题思路、代码3.1解题思路Gauss-Seidel迭代公式:x i(k+1)=(b i-∑-=1i i j a ij x j(k+1)-∑+=nij1a ij x j(k))/a ij(i=1,2,…,n)3.2 代码function x = GaussSeidel(A, b, es, maxit)% GaussSeidel: Gauss Seidel method% x = GaussSeidel(A, b):Gauss Seidel without relaxation% input:% A = coefficient matrix% b = right hand side vector% es = stop criterion(default = 0.00001%)% maxit = max iteration (default = 50)% output:% x = solution vectorif nargin < 2, error('at least 2 input arguments required'), end if nargin<4 | isempty(maxit), maxit=50; endif nargin<3 | isempty(es), es=0.00001; endk=0xk=[0 0 0 0][m, n] = size(A);if m~=n, error('Matrix A must be square'); endC = A;for i = 1:nC(i,i) = 0;x(i) = 0;endx = x';for i = 1:nC(i,1:n) = C(i,1:n)/A(i,i);endfor i = 1:nd(i) = b(i)/A(i,i);enditer = 0;while(1)xold = x;for i = 1:nx(i) = d(i)-C(i,:)*x;if x(i) ~= 0ea(i) = abs((x(i)-xold(i))/x(i)) * 100;endendk=k+1xk=x'%此行不打分号,并且转置,以便于输出每次迭代的结果 iter=iter + 1;if (max(ea)<=es | iter == maxit) break; end endend4、实验步骤4.1输入:4.2输出:……………….5、讨论和分析GaussSeidel迭代法是通过利用x i(k+1)=(b i-∑-=1i i j a ij x j(k+1)-∑+=nij1a ij x j(k))/a ij(i=1,2,…,n)这个公式,经过若干次运算,使结果越来越逼近方程的真实解。
迭代法求解方程的MATLAB实现1.引言迭代法是一种求解方程的常用方法,尤其适用于大规模矩阵和高维问题。
在迭代法中,我们通过不断迭代来逐步逼近方程的解。
本篇文章将介绍如何使用MATLAB实现迭代法求解方程。
2.收敛性判断在使用迭代法求解方程时,我们需要判断迭代是否收敛。
通常,我们使用以下两种方法进行收敛性判断:2.1 判断迭代公式是否收敛对于许多迭代公式,我们可以根据其结构来判断其是否收敛。
例如,Jacobi迭代法和Gauss-Seidel方法通常适用于对角占优的矩阵,而SOR方法适用于对角占优或松弛型的矩阵。
2.2 判断迭代误差是否收敛我们还可以通过判断迭代误差是否收敛来判断迭代是否收敛。
迭代误差通常定义为实际解与迭代解之间的范数。
如果迭代误差逐渐减小并趋于零,则说明迭代收敛。
3.迭代公式下面我们以Jacobi迭代法为例,介绍迭代公式的实现。
Jacobi迭代法的迭代公式如下:x{n+1}=(\frac{1}{a{ii}})(b i-A{ii,1:i-1}x1-A{ii,i+1:n}x_n)其中,A是系数矩阵,b是右侧向量,x是解向量,a_{ii}是矩阵A的主对角线元素。
4.误差计算为了判断迭代是否收敛,我们需要计算迭代误差。
通常,我们使用实际解与迭代解之间的相对误差或范数误差来衡量误差大小。
例如,相对误差可以按下式计算:||x^-x_n||_2/(||x^||_2)其中,x^*是实际解向量,x_n是第n次迭代的解向量。
5.MATLAB实现下面是一个使用MATLAB实现Jacobi迭代法的示例代码:function x = jacobi(A, b, x0, tol, max_iter)% 输入参数:系数矩阵A、右侧向量b、初始解向量x0、容许误差tol和最大迭代次数max_iter% 输出参数:方程的解向量xn = length(b); % 方程的未知数个数x = zeros(n, 1); % 初始化解向量xx(:) = x0; % 将初始解向量赋值给xerr = tol + 1; % 初始化误差大于容许误差,表示未收敛k = 0; % 初始化迭代次数k=0while err > tol && k < max_iterk = k + 1; % 更新迭代次数k=k+1for i = 1:n % 对每个未知数进行更新x(i) = (b(i) - A(i, 1:i-1)*x(1:i-1) - A(i, i+1:n)*x(i+1:n)) / A(i, i); % 更新解向量x的第i个元素x(i)的公式为x(i)=[b(i)-A(i,1:i-1)*x(1:i-1)-A(i,i+1:n)*x(i+1:n)]/A(i,i)endforendwhileif k < max_itererr = norm(x - x_prev, 'fro') / norm(x_prev, 'fro'); % 计算相对误差endendendfunction```。
解(1):采用Jacobi迭代法时,Matlab计算程序为: clear clci=1;a=[5 2 1;-1 4 2;2 -3 10];d=diag(diag(a));l=d-tril(a);u=d-triu(a);d0=inv(d);b=[-12;20;3];x0=[1;1;1];B=d0*(l+u);f=d0*b;x=B*x0+f;while norm(x-x0,inf)>=1e-4x0=x;x=B*x0+f;i=i+1;endxi采用Gauss-Seidel迭代法计算时,Matlab计算程序为: clearclci=1;a=[5 2 1;-1 4 2;2 -3 10];d=diag(diag(a));l=d-tril(a);u=d-triu(a);b=[-12;20;3];x0=zeros(3,1);B=inv(d-l)*u;f=inv(d-l)*b;x=B*x0+f;while norm(x-x0,inf)>=1e-4x0=x;x=B*x0+f;i=i+1;endxi习题6.7function [n,x]=sor22(A,b,X,x1,nm,w,ww)%用超松弛迭代法求解方程组Ax=b%输入:A为方程组的系数矩阵,b为方程组右端的列向量,X为迭代初值构成的列向量,x1为方程的精确解,nm为最大迭代次数,w为误差精度,ww为松弛因子%输出:x为求得的方程组的解构成的列向量,n为迭代次数n=1;m=length(A);D=diag(diag(A)); %令A=D-L-U,计算矩阵DL=tril(-A)+D; %令A=D-L-U,计算矩阵LU=triu(-A)+D; %令A=D-L-U,计算矩阵UM=inv(D-ww*L)*((1-ww)*D+ww*U); %计算迭代矩阵g=ww*inv(D-ww*L)*b; %计算迭代格式中的常数项%下面是迭代过程while n<=nmx=M*X+g; %用迭代格式进行迭代if norm(x1-X,'inf')<wdisp('迭代次数为');ndisp('方程组的解为');xreturn;%上面:达到精度要求就结束程序,输出迭代次数和方程组的解endX=x;n=n+1;end%下面:如果达到最大迭代次数仍不收敛,输出警告语句及迭代的最终结果(并不是方程组的解)disp('在最大迭代次数内不收敛!');disp('最大迭代次数后的结果为');xa=[4 -1 0;-1 4 -1;0 -1 4];b=[1;4;-3];c=200;d=5e-3;f=1.03;k=[0 ;0; 0];x1=[1/2;1;-1/2];g=sor22(a,b,k,x1,c,d,f)习题6.8function [n,x]=sor(A,b,X,nm,w,ww)%用超松弛迭代法求解方程组Ax=b%输入:A为方程组的系数矩阵,b为方程组右端的列向量,X为迭代初值构成的列向量,nm为最大迭代次数,w为误差精度,ww为松弛因子%输出:x为求得的方程组的解构成的列向量,n为迭代次数n=1;m=length(A);D=diag(diag(A)); %令A=D-L-U,计算矩阵DL=tril(-A)+D; %令A=D-L-U,计算矩阵LU=triu(-A)+D; %令A=D-L-U,计算矩阵UM=inv(D-ww*L)*((1-ww)*D+ww*U); %计算迭代矩阵g=ww*inv(D-ww*L)*b; %计算迭代格式中的常数项%下面是迭代过程while n<=nmx=M*X+g; %用迭代格式进行迭代if norm(x-X,'inf')<wdisp('迭代次数为');ndisp('方程组的解为');xreturn;%上面:达到精度要求就结束程序,输出迭代次数和方程组的解endX=x;n=n+1;end%下面:如果达到最大迭代次数仍不收敛,输出警告语句及迭代的最终结果(并不是方程组的解)disp('在最大迭代次数内不收敛!');disp('最大迭代次数后的结果为');xa=[5 2 1;-1 4 2;2 -3 10];b=[-12;20;3];c=200;d=5e-6;f=0.9;k=[0;0;0];g=sor(a,b,k,c,d,f)。
解(1):
采用Jacobi迭代法时,Matlab计算程序为:clear
clc
i=1;
a=[5 2 1;-1 4 2;2 -3 10];
d=diag(diag(a));
l=d-tril(a);
u=d-triu(a);
d0=inv(d);
b=[-12;20;3];
x0=[1;1;1];
B=d0*(l+u);
f=d0*b;
x=B*x0+f;
while norm(x-x0,inf)>=1e-4
x0=x;
x=B*x0+f;
i=i+1;
end
x
i
采用Gauss-Seidel迭代法计算时,Matlab计算程序为:
clear
clc
i=1;
a=[5 2 1;-1 4 2;2 -3 10];
d=diag(diag(a));
l=d-tril(a);
u=d-triu(a);
b=[-12;20;3];
x0=zeros(3,1);
B=inv(d-l)*u;
f=inv(d-l)*b;
x=B*x0+f;
while norm(x-x0,inf)>=1e-4
x0=x;
x=B*x0+f;
i=i+1;
end
x
i
习题6.7
function [n,x]=sor22(A,b,X,x1,nm,w,ww)
%用超松弛迭代法求解方程组Ax=b
%输入:A为方程组的系数矩阵,b为方程组右端的列向量,X为迭代初值构成的列向量,x1为方程的精确解,nm为最大迭代次数,w为误差精度,ww为松弛因子
%输出:x为求得的方程组的解构成的列向量,n为迭代次数
n=1;
m=length(A);
D=diag(diag(A)); %令A=D-L-U,计算矩阵D
L=tril(-A)+D; %令A=D-L-U,计算矩阵L
U=triu(-A)+D; %令A=D-L-U,计算矩阵U
M=inv(D-ww*L)*((1-ww)*D+ww*U); %计算迭代矩阵
g=ww*inv(D-ww*L)*b; %计算迭代格式中的常数项
%下面是迭代过程
while n<=nm
x=M*X+g; %用迭代格式进行迭代
if norm(x1-X,'inf')<w
disp('迭代次数为');n
disp('方程组的解为');x
return;
%上面:达到精度要求就结束程序,输出迭代次数和方程组的解
end
X=x;n=n+1;
end
%下面:如果达到最大迭代次数仍不收敛,输出警告语句及迭代的最终结果(并
不是方程组的解)
disp('在最大迭代次数内不收敛!');
disp('最大迭代次数后的结果为');x
a=[4 -1 0;-1 4 -1;0 -1 4];
b=[1;4;-3];
c=200;
d=5e-3;
f=1.03;
k=[0 ;0; 0];
x1=[1/2;1;-1/2];
g=sor22(a,b,k,x1,c,d,f)
习题6.8
function [n,x]=sor(A,b,X,nm,w,ww)
%用超松弛迭代法求解方程组Ax=b
%输入:A为方程组的系数矩阵,b为方程组右端的列向量,X为迭代初值构成的列向量,nm为最大迭代次数,w为误差精度,ww为松弛因子
%输出:x为求得的方程组的解构成的列向量,n为迭代次数
n=1;
m=length(A);
D=diag(diag(A)); %令A=D-L-U,计算矩阵D
L=tril(-A)+D; %令A=D-L-U,计算矩阵L
U=triu(-A)+D; %令A=D-L-U,计算矩阵U
M=inv(D-ww*L)*((1-ww)*D+ww*U); %计算迭代矩阵
g=ww*inv(D-ww*L)*b; %计算迭代格式中的常数项
%下面是迭代过程
while n<=nm
x=M*X+g; %用迭代格式进行迭代
if norm(x-X,'inf')<w
disp('迭代次数为');n
disp('方程组的解为');x
return;
%上面:达到精度要求就结束程序,输出迭代次数和方程组的解
end
X=x;n=n+1;
end
%下面:如果达到最大迭代次数仍不收敛,输出警告语句及迭代的最终结果(并不是方程组的解)
disp('在最大迭代次数内不收敛!');
disp('最大迭代次数后的结果为');x
a=[5 2 1;-1 4 2;2 -3 10];
b=[-12;20;3];
c=200;
d=5e-6;
f=0.9;
k=[0;0;0];
g=sor(a,b,k,c,d,f)。