近似解 牛顿法 初值x=1 1.9449e-04 迭代次数n 21 不动点法 0.0182 k=410
(3)为了提高牛顿法求重根的收敛阶采用以下俩种方法： 方法一：
程序：function[p1,err,y]=newtonroot(f,df,p0,eps,max1) p0
for k=1:max1
p1=p0-2*（feval('f',p0)/feval('df',p0)）； err=abs(p1-p0); p0=p1;
p1,err,k,y=feval('f',p1) if(err end 取初值p=0.5,当迭代次数k=9满足精度要求,近似解为2.4925e-05 方法二：程序： function[p1,err,y]=newtonroot(f,df,ddf,p0,eps,max1) p0 for k=1:max1 p1=p0-(feval('f',p0)*feval('df',p0))/((feval('df',p0))^2-feval('f',p0)*feval('ddf',p0)) err=abs(p1-p0); p0=p1; p1,err,k,y=feval('f',p1) if(err 取初值p=0.5,k =3执行结果为： p1 =1.5022e-08 err =2.2741e-08 ans = 1.5022e-08 以上四种方法中简单迭代法的收敛速度最慢 最后一种收敛速度最快。 3-4 试验目的：体验Steffensen’s method加速技巧 试验内容： 先用Newton法求解方程x?tanx?0，再用Steffensen’s method求解，比较迭代步数，精确到0.00001。 （1）简单原理：其基本思想是：对给定的初值p0(0),首先应用不动点迭代法pn+1=g(pn)计算俩步，然后用Atiken公式加速。 （2）程序代码 Newton法迭代 程序：function[p1,err,y]=newtonroot(f,df,p0,eps,max1) p0 for k=1:max1 p1=p0-feval('f',p0)/feval('df',p0); err=abs(p1-p0); p0=p1; p1,err,k,y=feval('f',p1) if(err 定义迭代函数：存储为f.m function y=f(x) y=x-tan(x) end 定义迭代函数：存储为df.m function g=df(x) g=1-1/(cos(x))^2 end Atiken公式加速程序： function [n,err,p1,p2]=steffen(g,p0,tol,max1) p0 n=0,p(1)=p0; while n<=max1 for k=2:3 p(k)=feval('g',p(k-1)); end p1=p(1)-(p(2)-p(1))^2/(p(3)-2*p(2)+p(1)); err=abs(p1-p(1)); n=n+1; p(1)=p1 n,err; if(err （3）结果与分析 牛顿法 取初值p0=0.5，迭代次数k=20满足精度要求. p1 =1.6018e-04 err =8.0089e-05 k =20 Ans=1.6018e-04 Atiken公式加速 取初值p0=0.5，迭代次数k=20满足精度要求。 err =8.6054e-05 n =20 p =0.1776×1.0e-03 二者迭代步数相同。从结果上看牛顿法更精确一些。 上机作业2 数值试验5-1 实验目的：熟悉Jacobi、Seidel、Sor迭代法，了解松弛因子对收敛速度的影响。 实验内容：分别用Jacobi、Seidel、Sor???0.8,1.1,1.2,1.3,1.4,1.5?迭代法求解下面的方程组，并作结果分析。 初值x(0)??0,0,0,0,0?，精度要求：x?k?1??x(k)T??10?5x(k?1)? （1） （2） ?12.3x1?2x2?x3?3.4x4?3.7x5?4.8?13.3x1?4x2?x3?3.5x4?3.8x5?5.8?1.4x?9x?3x?2.4x?2.7x?2.3?3.4x?9x?3x?4.4x?2.3x?4.31234512345?????2.1x1?x2?8x3?2.6x4?5.8x5?2.5?4.1x1?x2?7x3?2.7x4?5.9x5?2.6?3.5x?2.1x?x?13x?4.6x?3.6?2.5x?2.4x?x?13x?5.6x?3.81234512345????2.5x1?x2?2x3?5.3x4?14.8x5?2.2??1.5x1?x2?3x3?4.3x4?14.9x5?4.2 1(1)实验原理 Jacobi迭代法求解方程组（1） 程序代码： 先建立M文件，保存文件名为jacobi.m。 function y=jacobi(a,b,x0) D=diag(diag(a)); U=-triu(a,1); L=-tril(a,-1); B=D\\(L+U); f=D\\b; y=B*x0+f;n=1; while norm(y-x0)>=1.0e-5 x0=y; y=B*x0+f;n=n+1; end y n Seidel迭代法 D=diag(diag(A)); D0=inv(D); L=tril(A)-D; U=triu(A)-D; M=L+D; M0=inv(M); B=-M0*U; d=M0*b; x0=[0 0 0 0 0]'; x=B*x0+d; while norm(x-x0)>1e-5*norm(x) x0=x; x=B*x0+d; i=i+1; end Sor迭代式求解方程组 程序代码： 建立文件名为sor.m的M文件，具体如下： function y=sor(a,b,w,x0) D=diag(diag(a)); U=-triu(a,1); L=-tril(a,-1); B=(D-w*L)\\((1-w)*D+w*U); f=(D-w*L)\\b*w; y=B*x0+f;n=1; while norm(y-x0)>=1.0e-5 x0=y; y=B*x0+f; n=n+1; End y n 2程序计算结果 三种迭代法求解方程组（1）结果列表 Sor 0.8 1.1 1.2 1.3 1.4 1.5 NaN NaN NaN Jacobi — 0.3906 0.1678 0.0996 Seidel — 0.3906 0.1678 0.0996 ?x1?0.3906 0.3906 0.3906 0.3906 0.3906 ???x2?0.1678 0.1678 0.1678 0.1678 0.1678 x*??x3?0.0996 0.0996 0.0996 0.0996 0.0996 ?? ?x4??x??5?
